{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# _*Pricing Iron Condor Option*_ \n",
    "\n",
    "The latest version of this notebook is available on https://github.com/Qiskit/qiskit-iqx-tutorials.\n",
    "\n",
    "***\n",
    "### Contributors\n",
    "Stefan Woerner<sup>[1]</sup>, Daniel Egger<sup>[1]</sup>\n",
    "### Affliation\n",
    "- <sup>[1]</sup>IBMQ"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Introduction\n",
    "<br>\n",
    "Suppose a <a href=\"http://www.theoptionsguide.com/iron-condor.aspx\">iron condor option</a> with strike prices $K_1 < K_2 < K_3 < K_4$ and an underlying asset whose spot price at maturity $S_T$ follows a given random distribution.\n",
    "The corresponding payoff function is defined as:\n",
    "<br>\n",
    "<br>\n",
    "$$ F(S_T) = \n",
    "\\begin{cases}\n",
    "0         ,& S_T < K_1 \\\\\n",
    "S_T - K_1 ,& K_1 \\leq S_T < K_2 \\\\\n",
    "K_2 - K_1 ,& K_2 \\leq S_T < K_3 \\\\\n",
    "K_3 - S_T ,& K_3 \\leq S_T < K_4 \\\\\n",
    "0         ,& S_T \\geq K_4. \n",
    "\\end{cases}$$\n",
    "<br>\n",
    "In the following, a quantum algorithm based on amplitude estimation is used to estimate the expected payoff, i.e., the fair price before discounting, for the option:\n",
    "<br>\n",
    "<br>\n",
    "$$\\mathbb{E}\\left[ F(S_T) \\right].$$\n",
    "<br>\n",
    "The approximation of the objective function is explained in detail in the following paper:<br>\n",
    "<a href=\"https://arxiv.org/abs/1806.06893\">Quantum Risk Analysis. Woerner, Egger. 2018.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import numpy as np\n",
    "\n",
    "from qiskit import BasicAer\n",
    "from qiskit.aqua.algorithms import AmplitudeEstimation\n",
    "from qiskit.aqua.components.uncertainty_models import LogNormalDistribution\n",
    "from qiskit.aqua.components.uncertainty_problems import UnivariateProblem\n",
    "from qiskit.aqua.components.uncertainty_problems import UnivariatePiecewiseLinearObjective as PwlObjective"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Uncertainty Model\n",
    "\n",
    "We construct a circuit factory to load a log-normal random distribution into a quantum state.\n",
    "The distribution is truncated to a given interval $[low, high]$ and discretized using $2^n$ grid points, where $n$ denotes the number of qubits used.\n",
    "The unitary operator corresponding to the circuit factory implements the following: \n",
    "$$\\big|0\\rangle_{n} \\mapsto \\big|\\psi\\rangle_{n} = \\sum_{i=0}^{2^n-1} \\sqrt{p_i}\\big|i\\rangle_{n},$$\n",
    "where $p_i$ denote the probabilities corresponding to the truncated and discretized distribution and where $i$ is mapped to the right interval using the affine map:\n",
    "$$ \\{0, \\ldots, 2^n-1\\} \\ni i \\mapsto \\frac{high - low}{2^n - 1} * i + low \\in [low, high].$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# number of qubits to represent the uncertainty\n",
    "num_uncertainty_qubits = 3\n",
    "\n",
    "# parameters for considered random distribution\n",
    "S = 2.0 # initial spot price\n",
    "vol = 0.4 # volatility of 40%\n",
    "r = 0.05 # annual interest rate of 4%\n",
    "T = 40 / 365 # 40 days to maturity\n",
    "\n",
    "# resulting parameters for log-normal distribution\n",
    "mu = ((r - 0.5 * vol**2) * T + np.log(S))\n",
    "sigma = vol * np.sqrt(T)\n",
    "mean = np.exp(mu + sigma**2/2)\n",
    "variance = (np.exp(sigma**2) - 1) * np.exp(2*mu + sigma**2)\n",
    "stddev = np.sqrt(variance)\n",
    "\n",
    "# lowest and highest value considered for the spot price; in between, an equidistant discretization is considered.\n",
    "low  = np.maximum(0, mean - 3*stddev)\n",
    "high = mean + 3*stddev\n",
    "\n",
    "# construct circuit factory for uncertainty model\n",
    "uncertainty_model = LogNormalDistribution(num_uncertainty_qubits, mu=mu, sigma=sigma, low=low, high=high)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZwAAAEyCAYAAADOV2anAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3debgcVZnH8e9PEAgEMLKERSSACiJRMEGMoCSArDpBEKKiYxCNjAKOIgwiQgAdAWVxYBQjakRGoyKiIIhhSRDZ40JYomwhBhQEAjGEQELe+eOcC51K973dt/tWNcnv8zz93FtVp6re7tu33646myICMzOzgfaKqgMwM7OVgxOOmZmVwgnHzMxK4YRjZmalcMIxM7NSOOGYmVkpKk84kraVdI2khZIekXSKpFX62OdNkn6Tyz8naY6kCyRtXCg3WVLUeWwzsM/KzMyKVq3y5JKGAFcDdwNjga2AM0mJ8IRedl0XeBC4EHgE2AI4CRghaceIWFJTdhZwaGH/2c3Et/7668ewYcOaKTognnnmGdZaa63Kzt9It8YFjq0/ujUucGz9UXVcM2bMeDwiNqi7MSIqewBfAOYB69SsOxZYWLuuyWO9GwjgrTXrJgO39ze+ESNGRJWuu+66Ss/fSLfGFeHY+qNb44pwbP1RdVy9feZWfUttH+CqiJhfs24KMAjYtcVjPZF/rtaJwMzMrLOqTjjbkG55vSgi5pCucPqsZ5H0CkmrSdoaOA24Dbi1UGxbSfNzXc8NklpNZGZm1gGKCsdSk7QYOCYizimsnwtcGBHH97H/b4C98uIMYN+IeKxm+2eA50l1RBsARwMjgF0iopiYevaZAEwAGDp06IgpU6b056l1xIIFCxg8eHBl52+kW+MCx9Yf3RoXOLb+qDquMWPGzIiIkXU3NrrXVsYDWAz8Z531c4H/bmL/1wM7AR8mXSnNANbopfyapMYGlzYTn+tw6uvWuCIcW390a1wRjq0/qo6LLq7DmUdqcVY0JG/rVUTcGxG3RMRFpCudHYAP9VJ+IXAF8Nb+hWtmZv1VdcKZRaGuRtJmpCuRWXX3aCAiHgKeBLbsq2h+mJlZiapOOFcCe0lau2bdOOBZYHorB8oNB9Yj3TJrVGYQsB/p1puZmZWo0o6fwPnAUcAlkk4nXZ1MBM6KmqbSku4DpkfEYXn568AS4BbgKeCNpP4795OaVSNpXeBy4CLgPmB94LPAJsBBJTw3MzOrUWnCiYh5knYHzgMuIyWPs0lJp9aqQO1wN7cDR5Jak60BzAF+Dnw1Ip7JZZ4D/kkasWBDYBFwE7BrRNw+EM/HzMwaq/oKh4i4G9itjzLDCstTyFcyveyzCDig3fjMAIYd9+u2j3H08CWMb/M4s0/br+04zKpSdR2OmZmtJJxwzMysFE44ZmZWCiccMzMrhROOmZmVwgnHzMxK4YRjZmalcMIxM7NSOOGYmVkpnHDMzKwUTjhmZlYKJxwzMyuFE46ZmZXCCcfMzErhhGNmZqVwwjEzs1I44ZiZWSmccMzMrBSVJxxJ20q6RtJCSY9IOkXSKn3s8yZJv8nln5M0R9IFkjauU3aspJmSFkm6W9K4gXs2ZmbWyKpVnlzSEOBq4G5gLLAVcCYpEZ7Qy67rAg8CFwKPAFsAJwEjJO0YEUvy8XcBfg58EzgK2Bf4saR5EfHbAXlSZmZWV6UJBzgcGAQcEBHzgamS1gEmSjojr1tORNwI3FizapqkucBvgTcDf8jrvwRcHxFH5eXrJL0JODGXNTOzklR9S20f4KpCYplCSkK7tnisJ/LP1QAkrQ6MAX5aKDcFGCVp3dbDNTOz/qo64WwDzKpdERFzgIV5W68kvULSapK2Bk4DbgNuzZu3Al5ZPD5wD+l5v6G90M3MrBWKiOpOLi0GjomIcwrr5wIXRsTxfez/G2CvvDgD2DciHsvbdgZuAHaIiD/V7PM64F5gr3r1OJImABMAhg4dOmLKlCn9fXptW7BgAYMHD67s/I10a1wwcLHNfPjpto8xdBA8+mx7xxi+aecvzFfGv2cndGtsVcc1ZsyYGRExst62qutw2nUk8Grg9aRGBldK2jkiFvX3gBExCZgEMHLkyBg9enQn4uyXadOmUeX5G+nWuGDgYht/3K/bPsbRw5dw5sz2/uVmHzK67TiKVsa/Zyd0a2zdGhdUn3DmkVqcFQ3J23oVEffmX2+R9DtSy7UPAd+r2b94/CE15zYzs5JUXYczi0JdjaTNgDVZvu6lVxHxEPAksGVedT+wuHj8vLwU+Gs/4jUzs36qOuFcCewlae2adeOAZ4HprRwoNxxYj3SVQ0Q8B1wHHFQoOg64KSLavylvZmZNq/qW2vmkDpmXSDqddHUyETirtqm0pPuA6RFxWF7+OrAEuAV4CngjcCzpqqa2lv9UUh+dc4BLSR0/9wX2HtinZWZmRZUmnIiYJ2l34DzgMlLyOJuUdGqtCtQOd3M7qcHABGANYA5pRIGvRsQzNce/QdL7gS8D/0Gu4/EoA7aiGNahxgztNoqYfdp+bcdhK76qr3CIiLuB3fooM6ywPIVlr2R62/dS0tWNmZlVqOo6HDMzW0k44ZiZWSmccMzMrBROOGZmVgonHDMzK4UTjpmZlcIJx8zMSuGEY2ZmpXDCMTOzUjjhmJlZKZxwzMysFE44ZmZWCiccMzMrRcujRUsaDrwN2Ig0NcCTpNkzb4wIT9tsZmZ1NZVwJG1Jmk/mEGAoaYrmp4DngFeRpoReKmk6cAHwk4hYOiARm5nZy1Kft9QkXQDcBWwPnALsAKwRERtExGsiYjCwIfBeYCZwBnCPpF0GLmwzM3u5aeYK51lgm4h4qFGBiHgcuBK4UtLngIOATTsTopmZrQj6vMKJiCN7SzZ1yi+NiJ9ExE+aKS9pW0nXSFoo6RFJp0hapY99dpT0fUn35f3+IukkSWsUyk2UFHUeezf7fMzMrDPammJa0nbAroCA6RExs8X9hwBXA3cDY4GtgDNJifCEXnYdl8ueDtwLvBk4Nf88sFD2aaCYYO5pJU4zM2tfvxOOpP8AvgJcA6wFfE3S0RHxzRYOczgwCDggIuYDUyWtA0yUdEZeV89p+TZej2mSFgHflrR54YpsSUTc3EJMZmY2AJppNLBmg03/BYyKiIMiYl/gCOCLLZ5/H+CqQmKZQkpCuzbaqZBsevwx/9ykxRjMzKwEzXT8/KukQ+qsF6l5dI/+NIPeBphVuyIi5gAL87ZWjMox3F9Y/ypJj0taLOmPkg7oR5xmZtYmRUTvBaR3AecAzwNHRcStef2nSc2kryH1w9kdODYizm365NJi4JiIOKewfi5wYUQc3+RxNgLuAK6IiPE16z9MarL9R2Bt4JPAvsCBEXFJg2NNACYADB06dMSUKVOafTodt2DBAgYPHlzZ+Rvp1rhg4GKb+fDTbR9j6CB49Nn2jjF803WXWe7WuDplZXyvtavquMaMGTMjIkbW29ZnwgGQJOAwUsX8VOC/IuLvkt7CS7e+ro+IP7USWCcSjqTVSA0PXgOM6G20g/w8bgQGRcT2fR175MiRcfvtt/dVbMBMmzaN0aNHV3b+Rro1Lhi42IYd9+u2j3H08CWcObOtdjrMPm2/ZZa7Na5OWRnfa+2qOi5JDRNOU2OpRXIBsDXwKHCnpC8CsyLif/KjpWSTzQPqfTUakrf1KieQC4E3Afv2NbROpOx6CfDmvppem5lZZ7U0eGdEzI+IY4CdSOOpzZL0/jbOP4tCXY2kzUi36GbV3WNZ55CaU4+NiGbKA0R+mJlZiZpqpSbpy5JuyZXuk4BFETGWVNdxkqTp+fZaq64E9pK0ds26caTRDab3EdcXSC3jPhwRNzRzsnxFdCDw54h4oR/xmplZPzVz4/a7wLakPjcLSUlmqqRtI2JqTjSfyusujYgJLZz/fOAo4BJJpwNbAhOBs2qbSku6j9Sx9LC8/CHgv4HJwMOS3l5zzPsj4p+53HTg56SrpbWAT5CuzvZvIUYzM+uAZhLOPsBBETEVQNLvgSdIPf3vy6NCnyfpR8BJrZw8IuZJ2h04D7iMNAL12aSkU4yzts5lz/xzfH7UOpSUiADuA/4T2JjUZPoPwH4RcWUrcZqZWfuaSTizgI9ImgEsIjUtfgaYW1soIp4EPtNqABFxN7BbH2WGFZbHs3yiqbffYa3GY2ZmA6OZhPNR0hXD46TK9tmkK55FAxeWmZmtaPpMOBHxF2CUpLWA1Tyrp5mZ9UfTvb0i4hnSrTQzM7OWNdMs+iOtdpKU9DpJ7+x/WGZmtqJppuPn54D7JZ3aW18bSetJOkTSZcCfSC3DzMzMgObqcHaQNA44EviipAWkCcweB54DXgVsAbyWNBzNRcDhEfHwgEVtZmYvO03V4eTpon8iaStgD+CtwEakzpSPAtcDvwemRcTiAYrVzMxexloaIjYi7mf5+WbMzMz61NLgnWZmZv3lhGNmZqVwwjEzs1I44ZiZWSlaSjiS3ivJScrMzFrWavK4FJgr6XRJbxyIgMzMbMXUasLZCvgOcDBwp6SbJH1C0jqdD83MzFYkLSWciJgdESdFxBbAu0kTnJ0N/F3SDyWNGYggzczs5a/f9TERcW1EfAR4AzADOAS4WtIDkj4rqaVOpWZmtmLrd8KRtKukycBfgO2A/yVN/XwxcDJwYScCNDOzFUOrrdQ2l3SipPuBa4HNgAnAxhFxZERcExHHkmYJHdvkMbeVdI2khZIekXRKX9MhSNpR0vcl3Zf3+4ukkyStUafszpJukbRI0oOSjmrlOZuZWWe0etvrAeAR0pTT34uIBxuUuwu4ta+DSRoCXA3cTUpQWwFnkhLhCb3sOi6XPR24F3gzcGr+eWDN8V8HXAVcDnwBeBtwlqSFEXFBX/GZmVnntJpw3gNcFRFLeysUEX8FmmlAcDgwCDggIuYDU3OLt4mSzsjr6jktIh6vWZ4maRHwbUmbR8RDef0xpAT54YhYAlwr6bXASZK+GxHRRIxmZtYBrdbh7EialmA5kjaWdGKLx9uHlMBqE8sUUhLatdFOhWTT44/55yaF41+Sk03t8V9DqncyM7OStJpwTiJ9WNezSd7eim2AWbUrImIOsDBva8UoYCl5+gRJa5HqmGYVyt1Tc24zMyuJWrmrJGkpsFNE3FZn21jguxGxfgvHWwwcExHnFNbPBS6MiOObPM5GwB3AFRExPq/bFJgLvC8iLq0puyqwGPhkREyqc6wJpIYQDB06dMSUKVOafTodt2DBAgYPHlzZ+Rvp1rhg4GKb+fDTbR9j6CB49Nn2jjF803WXWe7WuDplZXyvtavquMaMGTMjIkbW29ZnHY6kj5JanQEE8C1JxbqVNYDhwG/bCbQ/JK0G/BRYAHy23ePlJDQJYOTIkTF69Oh2D9lv06ZNo8rzN9KtccHAxTb+uF+3fYyjhy/hzJntdU+bfcjoZZa7Na5OWRnfa+3q1riguUYDC4En8u8CngaeLJR5HrgS+GaL558H1PtqNCRv65Ukkfr7vAnYOSJq93kq/ywef0jNuc3MrCR9JpyI+BnwMwBJ3wdO6aU5dKtmUahLkbQZsCbL173Ucw6pOfW7I6JYF/SMpL8Vj1+z3MzxzcysQ1odS+3QDiYbSFdFe0lau2bdOOBZYHpvO0r6AnAEqcnzDb0c/32FjqTjgL8Bd/Y7ajMza1nVc9ucDzwHXCJpj1xhPxE4q7apdB5R4Ls1yx8C/pt0O+1hSW+veWxQc/yvkVrV/VDSGEnHAp8kXaW5D46ZWYmaaTRwKzA+Iu6WdBup4UBDEfG2Zk8eEfMk7Q6cB1xGqnc5m5R0inHWXqXsmX+Oz49ah5JGQiAi7pO0N3AW6WrnH8DRHmXAzKx8zTQauIt0i6vn945eGUTE3cBufZQZVlgez/KJptG+N5CGtDEzswo102jg0Jrfxw9oNGZmtsKqug7HzMxWEs3U4fRZb1OrlTocMzNbeTRbh+MWXWZm1pZm6nDGlxCHmZmt4FyHY2Zmpai0H46Zma08Ku+HY2ZmKwf3wzEzs1K0PAlGnn9mPKn3/sbA34FbgB9ExPMdjc7MzFYYLTUakPRG4F7gf4HtgBfyz/8F7pO0bccjNDOzFUKrVziTSBOwvTMi5vSslPRa4HLS6M/v6lx4Zma2omg14YwEPlibbAAiYo6kk4AfdSwyW+kM69B0ye1Ouzz7tP3ajsPMltdqP5zZwBoNtq0BzGmwzczMVnKtJpzjgC9L2ql2paS3A6cC/9WpwMzMbMXSn8E71wFulPQY8BiwYX48ARwPXDoAcZqZ2ctcfwbvvGuAYjEzsxVY5YN35qbU5wKjSFNMXwCcHBEv9LLPasBXgLeTGjKsERGqU24y8NE6h3hjRMxqP3ozM2tWyx0/O0nSEOBq4G5gLLAVcCapbumEXnZdE/g4cCtwI71PUT0LOLSwbnb/IjYzs/6qNOEAhwODgAMiYj4wVdI6wERJZ+R1y4mIpyS9OiJC0hH0nnCeiYibOx+6mZm1ouXpCSSNk3S1pDmSHis+WjzcPsBVhcQyhZSEdu1tx4jwIKJmZi8jrQ5t8yHgB8B9wGuAX5FGGHgFMB84r8Xzb0O65fWi3Kl0Yd7WCdtKmi/pOUk3SOo1kZmZ2cBQKxcKkv4IXAycBiwGRkbEHyStDUwFLo6Ir7dwvMXAMRFxTmH9XODCiDi+iWMcAZzboNHAZ4DnSXVEGwBHAyOAXSLi1gbHmwBMABg6dOiIKVOmNPt0Om7BggUMHjy4svM3MlBxzXz46baPMXQQPPps3+V6M3zTdZdb162xdWtcndKt/wPQvbFVHdeYMWNmRMTIettarcN5PfD7iHhB0gukPjlExL8knQ6cDTSdcAZaRHyjdlnSFaRm3ccD+zfYZxJpzDhGjhwZo0ePHuAoG5s2bRpVnr+RgYqr3SFpIA1tc+bM9qomZx8yerl13Rpbt8bVKd36PwDdG1u3xgWt1+HMB1bPvz8MvLFmm4D1WjzePKDeV6MheVtHRcRC4ArgrZ0+tpmZ9a7VrzW3AW8GriLV35woaQnpttWJQKutwWZRqKuRtBmp2fNA9ZMJPGupmVnpWk04XwU2z7+fmH//FulK6Tbgky0e70rgGElrR8S/8rpxpCmtp7d4rD5JGgTsB8zo9LHNzKx3LSWc3J/l5vz7U8BYSasDqzfqM9OH84GjgEtyHdCWwETgrNrjSboPmB4Rh9Ws2wdYC9g+L78/b7otIh6StC6pBd1FpFZ16wOfBTYBDupHrGZm1oaOTTEtqeUppiNinqTdSc2pLyMNbXM2KekU41ylsO5bvHS1BfCz/PNQYDLwHPBP0ogFGwKLgJuAXSPi9lbiNDOz9rWUcPIU078hXSXMII0WvR3w78CXJO0dEXe3csxcvreRAoiIYc2sK2xfBBzQSixmZjZwPMW0mZmVotVm0SOBE+tNMQ2cBOzYqcDMzGzF4immzcysFK3eUjsOOFPSgxFxS8/KmimmP9/J4Mzs5WtYh0ZBaHc0hdmn7dd2HNYZnmLazMxK4SmmzcysFJVPMW1mZiuHfg0RK2kTYBTwatKttJsj4pFOBmZmZiuWVjt+rgKcC3yCZXv+vyBpEnBkRCztYHxmZraCaLVZ9MnAx0iNA4aRpoIelpc/xvJD0piZmQGt31L7d+CEwqyec4CvSQrSQJwndio4MzNbcbR6hbMhcEeDbXfk7WZmZstpNeH8FfhAg20fAP7SXjhmZraiavWW2peBKXmwzouBR0lXNQcBY2icjMzMbCXX6gRsP5X0FKnxwDeAVwKLSVMV7B0RUzsfopmZrQiaTjiSXkmadO3OiBgl6RWkWTQfd1NoMzPrSyt1OC8A1wLbAETE0oh4zMnGzMya0XTCyYnlXmCjgQvHzMxWVK22UvsicKKk4Z0KQNK2kq6RtFDSI5JOySMa9LbPapK+Jul3kp7NfYAalR0raaakRZLuljSuU7GbmVnzWm2ldgKwHvAnSQ+TWqkt82EfEW9r9mCShgBXA3cDY4GtgDNJifCEXnZdE/g4cCtwI7Bbg+PvAvwc+CapU+q+wI8lzYuI3zYbp5mZta/VhHMXcGcHz384aXicAyJiPjBV0jrAREln5HXLiYinJL06IkLSETRIOMCXgOsj4qi8fJ2kN5FGQ3DCMTMrUavNosd3+Pz7AFcVEssU4HRgV+CyXmJpeBsNQNLqpL5BRxU2TQG+L2ndiHi6X1GbmVnLmqrDkTRI0oGSjpb0IUlDO3T+bYBZtSsiYg6wMG9rx1akfkKzCuvvIT3vN7R5fDMza4H6uFBA0pakepZhNavnAwe3Ww8iaTFwTEScU1g/F7gwIo5v4hhHAOdGhArrdwZuAHaIiD/VrH8dqbXdXvXilzQBmAAwdOjQEVOmTGn9iXXIggULGDx4cGXnb2Sg4pr5cPsXnEMHwaPPtneM4Zuuu9y6bo2tW+OC7o6tE1a2/89mjRkzZkZEjKy3rZlbamcAS4F3kkYU2IJUCf/t/PsKJSImAZMARo4cGaNHj64slmnTplHl+RsZqLjGH/frto9x9PAlnDmzX/MKvmj2IaOXW9etsXVrXNDdsXXCyvb/2QnN3FIbRZqS4PcRsSgi7gE+CbxW0sZtnn8eUO/rx5C8rd1jU+f4QwrbzcysBM0knI2BBwrr7gdE+51AZ1Goq5G0GanZc7HupVX3k8Z5K9YFbUO6Yvtrm8c3M7MWNNvxs/eKnv67EthL0to168YBzwLT2zlwRDwHXEcaybrWOOAmt1AzMytXszdHr5K0pM76a4rrI6KVSdjOJzVbvkTS6cCWpGmqz6ptKi3pPmB6RBxWs24fYC1g+7z8/rzptoh4KP9+KjBN0jnApaSOn/sCe7cQo5mZdUAzCefkgTp5RMyTtDtwHqnPzVPA2aSkU2tVoDjczbeAzWuWf5Z/HgpMzse/ISeiLwP/ATwIfMijDJiZla/PhBMRA5Zw8vHvpvFIAT1lhjWzrsG+l5KubszMrEKtDt5pZmbWL044ZmZWCiccMzMrhROOmZmVwgnHzMxK4YRjZmalcMIxM7NSOOGYmVkpnHDMzKwUTjhmZlYKJxwzMyuFE46ZmZXCCcfMzErhhGNmZqVwwjEzs1I44ZiZWSmccMzMrBSVJxxJ20q6RtJCSY9IOkVScTrpevutK+n7kuZJelrS/0lar1BmsqSo89hm4J6RmZnV0+cU0wNJ0hDgauBuYCywFXAmKRGe0MfuPwXeAHwcWAqcTppK+p2FcrOAQwvrZrcTt5mZta7ShAMcDgwCDoiI+cBUSesAEyWdkdctR9IoYE9g14i4Pq97GLhF0h4RcXVN8Wci4uaBfRpmZtaXqm+p7QNcVUgsU0hJaNc+9nu0J9kARMStwIN5m5mZdZmqE842pFteL4qIOcDCvK3p/bJ76uy3raT5kp6TdIOk3hKZmZkNEEVEdSeXFgPHRMQ5hfVzgQsj4vgG+00l3Srbv7D+ImDLiHhHXv4M8DypjmgD4GhgBLBLviKqd+wJwASAoUOHjpgyZUobz7A9CxYsYPDgwZWdv5GBimvmw0+3fYyhg+DRZ9s7xvBN111uXbfG1q1xQXfH1gkr2/9ns8aMGTMjIkbW21Z1Hc6Aiohv1C5LugK4Czge2L/BPpOASQAjR46M0aNHD3CUjU2bNo0qz9/IQMU1/rhft32Mo4cv4cyZ7b2tZx8yerl13Rpbt8YF3R1bJ6xs/5+dUPUttXlAva8fQ/K2ju4XEQuBK4C3thCjmZl1QNUJZxaFOhdJmwFrUr+OpuF+WaO6nVqRH2ZmVqKqE86VwF6S1q5ZNw54Fpjex34bSdqlZ4WkkcCWeVtdkgYB+wEz2gnazMxaV3XCOR94DrhE0h65wn4icFZtU2lJ90n6bs9yRNwE/Ba4UNIBkvYH/g+4oacPTh6J4HeSPilpd0njgOuATYD/LusJmplZUmmjgYiYJ2l34DzgMuAp4GxS0qm1KlAc7mZcLvs9UuK8HDiqZvtzwD9JIxZsCCwCbiJ1Fr29o0/EzMz6VHkrtYi4G9itjzLD6qx7ijRkTXHYmp7ti4ADOhCima1ghnWoBV07LfFmn7Zf2zG83FR9S83MzFYSTjhmZlYKJxwzMyuFE46ZmZXCCcfMzErhhGNmZqVwwjEzs1I44ZiZWSmccMzMrBSVjzRg5eqGHtawcvayNlvZ+QrHzMxK4YRjZmalcMIxM7NSOOGYmVkpnHDMzKwUTjhmZlYKJxwzMyuFE46ZmZWi8o6fkrYFzgVGAU8BFwAnR8QLfey3LnAOsD8pcV4OHBURTxTKjQW+DLweeCAf+yedfh5mZu1a0TtmV3qFI2kIcDUQwFjgFOBo4OQmdv8pMBr4ODAe2BG4tHD8XYCfA9cB+wC/Bn4sac+OPAEzM2ta1Vc4hwODgAMiYj4wVdI6wERJZ+R1y5E0CtgT2DUirs/rHgZukbRHRFydi34JuD4ijsrL10l6E3Ai8NuBe1pmZlZUdR3OPsBVhcQyhZSEdu1jv0d7kg1ARNwKPJi3IWl1YAzpSqjWFGBUviVnZmYlqTrhbAPMql0REXOAhXlb0/tl99TstxXwyjrl7iE97zf0I14zM+snRUR1J5cWA8dExDmF9XOBCyPi+Ab7TQWeiYj9C+svAraMiHdI2hm4AdghIv5UU+Z1wL3AXhGx3G01SROACXlxa+Av/X6C7VsfeLzC8zfSrXGBY+uPbo0LHFt/VB3X5hGxQb0NVdfhdJ2ImARMqjoOAEm3R8TIquMo6ta4wLH1R7fGBY6tP7o1Lqj+lto8oF5dypC8rZ39en4Wyw0pbDczsxJUnXBmUairkbQZsCb162ga7pfV1u3cDyyuU24bYCnw137Ea2Zm/VR1wrkS2EvS2jXrxgHPAtP72G+j3M8GAEkjgS3zNiLiOVL/m4MK+44DboqIp9sPf8B1xa29Oro1LnBs/dGtcYFj649ujavyRgNDgLuBO4HTSQnjLOCciDihptx9wPSIOKxm3VWk0QM+T7piOR14LCLeWVNmF2AacB6pU+i+ufze9RoMmJnZwKn0Cici5gG7A6sAl5FGGDgbOKlQdNVcptY40lXQ94ALgRnA+wrHvwF4P7AHcBXwb8CHnGzMzMpX6RWOmZmtPKquwzEzs5WEE1I1XdsAABkPSURBVI6ZmZXCCcfMzErhkQa6gCSRGjzsB7wReDWp5d0/gJuByRFRSb+h3C9qX0DAzyLiCUmvIbX22wqYDUyKiJklxvRfwBVlnrNZkgYBq0bEv2rWbQAcAWxL+rv+Cfjmy6RpfiXy/8R7gbeSpi+5nfQ374pK5zyq/ePAbrlxUlUx7AasBvw6Ip7J77VPk1r8PkD633ykivjqcaOBiuU3yBXACFKCeR7YlPRPdiXpjbM1cGpEnFpybG8DpgJrAUuAJ4G9crwvAHcB2wEbAXtExO9Kimsp6fWZBfwImBIR95dx7r5IugK4NyI+k5dHkf6OS0ktKUX6Wz9P+rC6q6S4dgAGRcSNNev2Br7AS4nwz8DE2jIlxXYjcFhE3JOXh5CmDxkBLMjFBpO+fO1Vm8wHOK5P9bJ5EPA14BuksRmJiG+WERe8OCbkNcBmedWDpClbpgKvInV835rUp3FERMwtK7ZeRYQfFT6AH5PesMNr1m0C/Ab4eV7elfSP97GSY5tK6jz7KtLI2+cBc4FfAq/MZVYnfaBeV2JcS4GvkmZ5fY6U/G4DPgtsWvHf83FgbM3yzaQPhrVr1q1LatJ/VYlx3Qx8sWb5Y/l1vAb4InBC/lsvqY2/xL/n22qWv0v6crN3zbq9ScNRnV1yXC/kn/UetdteKPk1+ynpC8LrSHdEfpg/R27sea+RBvH8M/DtMmPrNe6qA1jZH6RptQ+ss35YfkNvnJePB/5ccmxPAPvULG+Y/7n2LJTbD3i8xLhe/IDK/2wT8gfnkvyYlte9uoK/50LgXTXLzxdfr5rX7JkS45pfGwdwH3BunXLnV/A+KyacfwL/Wafc54GHSozrF8DfgUPJd4Nqtr0qx/2usuIpnP8R4OCa5c1zPAcUyh0K/LWKGOs93Gigeq8gJZaiF0i3X3oGH72FaubwiTq/F+/DVnZfNiKejIhJEbE78BrSFOWrkT44/y6p/UniW3MnaeK/Ho+SkmLReqTkVJalheXNgYvrlLuYdCumSq8i1dkUzSDdvi1FRLwP+ChwDHBbnvLkxc1lxdHAENIt+B4P558PFco9QPq/6ApOONWbCnxZ0pY9K/I97P8hvaF6GgsMBsquZL4d+LykwZJeQbrKehj4D0mr5FhXBT5F+qCtVET8IyK+ERHvALYgjVixSclhnAYcJ+lj+bX5CvA1Se+WtJqk1XPdyVdZfjbagfQ74JCa5buAekPY78hLH15lOlDSp3K9yTyg3nwq65Ou1EoTaVSSN5NmCv61pCm50UzVHiN9aejxAvBt0hecWhsCpdR5NaXqS6yV/UH69nEXaWTr+0hjyz1LutVWezvrDOAnJcc2kvTPvzjH9ATwFlISfIA0HNGDpHqUMSXGtcwtmG57AB8nfTA+Ddyaf3+BdLvvhfz4BbBmiTENz3H8EHgbaSr2x0gJ8d2kCufTgEXUuZ1Vwt+z+PhenXLfBn5X4d91I9IwWguAM/PfsapbapfWe43qlDsXmFrVa1Z8uJVaF8hXCweTPszXICWeH0XEk5UGBuRvc+8hNaH/eUT8XdJGwLGkWy8PARdExB9KjOkk4DvRRc09iyStRxrv722kDyqRkvc9wOURMaOCmLYHvgXsRLolpLyp5/d5wCkR8Y2yY2uGpE8A90fEtRXHMYo05uPWwH5Rcqu+HMNQ0heWB/so9zlSndw15UTWOyccs5WMpDeSkk4xEd4YEYurjM1WbE44XUTSm0gTxNXOSjorSuqr0SpJr4iIYmV0ZSStQeqMuhS4r+oPz1yHsyU1HXkjYk6VMb3c5A6gRIUfVLkzryJiYc267ckdn6u4Wn25cqOBLpArmB8C7gB+RppAaVL+/Q5JsyUdWlFsB0i6VNIVkt6b142TNBtYLOmhfKujzJg+LOljNcurSjqN1HfjDlIDhiclHVdmXDXxjJD0K1Jl7T3A74GbgAclPSzpFElrVhFbN5K0Z2ESRiTtL+kPpPrD5yXdLmm/kuNaV9IvSHVf8yV9R9Iqkn4A/IH0/3mrpBskrV9mbM2SdKCkeq1gK+GEUzFJR5IqQy8HRpNalbwyPzYkdfq8HDhf0qdLju1gUjPZ9Un/+D/JyeWHpH4vR5E6mp0vaa8SQzue1OG0x+k5lq8C7yK9ZmcCJ0k6vsS4kLQn6TXZhHSf/1RSr/kXgImkCQYPBG7MrRHLjO09kq6RdI+kX0p6V50yO1XwAXUlaUinnhjeB1xCasBwXH4sBn6ZX9+ynAq8E/gcqaPsO0gtC3cjdUQdSqrf3CKXtT74llrFJD0AnB8RZ/RR7ljg8IjYsrdynSTpNmBGRByelw8hTXh3XkQcXVPu+8BmEbFHSXEtJLXgm56XHwO+UqzslvR54MiI2LzOYQYqthnAnRHx0cL6I0l9hLYk9RO6Ebg5InobPqWTcb2bNHrFzcAfgVHA9sA5wOd7bllJ2olUl1Oc8HAgY1sKvD0ibs3LfwAejoj3FspdAawVEbuWFNds0vvqO3l5B1JfoEMj4gc15T4BHB8RW5QRVz7n95osujkwusy/Z298hVO9jUlNZ/tyKyV2esu2ZtnOgZeTrryKnSkvIY3HVZanSVddPdYlDeFR9GfSVWKZtgUuqrP+IuC1wNYRsYj0Qf++OuUGyknAhRGxc0QcEREjgE8AnwQuyfVf3WI70lV/0STSYJ5l2YCX+sFBHjONNE5Zrfuo329oIH2U1JR9eB+P0r5sNcMJp3p/Bj6RO1bWlStOP0GqnyhTsOzU3j0DKT5VKLeA1Du8LL8idUhdLS9fDXywTrkPkvo4lekxUvP2oreQXs+ezrsP8dIoEmXYjkIijIjvkW4/vh24VlK9ERHKUnur5Wnqd1Z8hnI/sx4kvT493klq/PGOQrmdgbIbg9wLXBsRO/b2IN2O7BqenqB6R5Nuddwt6RLSCMg9H+jrklqtvY/UQXTvkmN7iPSN/SqAiHgh90G4p1BuK9KYU2X5Aqnn/J2SLiB1QD1d0nakcdREGl5mB9IQ92WaBJwqaTApET5P6r3/RdIApz19h7ak3A+pRaRRv5cRETPykC1XkW7zTSwxplpXSVqSf1+X9LebXiizDeW+z84HviFpOCkJHkx6730p/33/TLri+izl1+HczPKJr57a/laVc8KpWET8PjexPJY09MhmhSJ/I1Wqfi3KH4L/EgrjMEXELXXKfQAobU6QiHhS0ttJH+Kf46XbZqPy43nSkEHvjIjbyoorx/aVXCdxHHBiz2rSqOD/WVN0MfDfJYZ2B2l0gV8VN0TEAznpXAFMLjGmHifXWfdYnXUHkka0LkVEnJfvPHyQ1DDg2Ig4X9Jc0tBTPePhnQ98vay4snNJLeX6Mp1lx/arlBsNdJncXLbn9tRTtW3/u5Wk15JiLXWcq5rzD2PZToz3d0EfnFeSrvzWAB6o6rWpieeTpNZ9OzQawULSWqQhd/aICN9u70W+zb1+RPyz6lheTpxwzMysFL6l1iWUpnLehPTt/PE629cH9o2IC0sPro58D/sPwCFl37ZSl0/jrC6clrvb6WUyXXK+sqmd+noGKd7Sv7lLGkm6zSjSNPSzJL2FdIuy5312XkRcVXZsjfgKp2KSVie1Hjogr1pKGpH2c7UflhX1j9i3l81rAT8h1VXcCRARV5QUV1dO45xj6cppuZulNM7aQRFxSonn7MrpktWlU1/nWPYiNZZ5ktR6bwNgLKne9R5SX6sRpAYrB0bEpWXF1qsyhqT2o9fhw08ktUr7BGk6gKNIc1rcC7y+ptxOlD+NbVdOsUuXTuOcz9uV03K3EP+BFbzPunK6ZLp06ut83t+ThtZZJS8fn+P4bqHcD0kdjCt/b0V4iunKH6Rm0EcU1m0EXE+aandUXldFwrmdl6bY3bzweHP+hzy4Z12JcXXlNM75nN06Lfdrm3wcXsH7rCunS66TcLpi6ut8zqdJV8g9y0NyvLsVyu1JatBTWmy9PVyHU73NKHTojIh/SNqd9O3k6jykTJn9D3rsSLryOp3U7+XzkeffkNTTafEfEVGc1nag9UzjfH1e7pZpnHt047Tcs5s8p5os10kvl+mSu2Lq6+xZlu1X1fP7oEK5NUl9sLqCE071HgFez0sfngBEatb7AUnfIF06l95YINJXpEmSfgp8GZgp6dz8e5VOA/5P0t9Ir0vPNM5PkG6j9XT8LHsaZ3hpWu4bSMmudlruayN1nq1iWu5/AdcCF/RRbhdSn7AydfN0yQfmynnooqmvSbfUTpR0bz7310mzBf+XpOsj4l/5S+GxpITYHaq+xFrZH6TBMKf1UeYLlFxP0iCON5N68j8CfIZqp9jtummcc1zdOi33VNJQKH2Vq6IOpyunS6aLp74m1XfNrnmv30+6JdrzvzCTlJznAduXGVtvD1/hVO+bwDhJr44GHfIi4qtK8+W8u9zQlovjDmC0pA8AZ1DhkBkRcUGeq6RnGucn6YJpnCPi9jwUSnFa7nfx0rTcV1LytNykK+gJTZT7J4Wr7RJ8knTrpy8PkpJTKaL5zq+3k1psliYi7stDOe1MapxyTUQ8K2k06cvY1qRb8j+Kklr1NcPNoq1f8m2htYAFEdE1EzyZWfdywjEzs1J4vKSXiTy97XerjqOebo2tW+OC7o6tW0m6WtI1VcdR1K1xQffF5jqcl48xdO8XhG6NrVvjgi6NTdLVpDsfu1cdSx2iC18zujcu6LLYfEvNzF6Uvw2/IiK6Zkh7W3F0Teaz3klaI08D0HW6NbZujQu6N7aI2L1bk42kV3bja9atcUH3xeaE8/KxH6lZaDfq1ti6NS7o0tiq+oCS9GlJ90t6VtKfJX2kTrG3UvJr1q1xdXtsjTjhmK0kuvUDKvfrOpc0COuXSJ0YJ0u6WNIaZcbycoir22PrjRsNVEzStU0WrTekxoDq1ti6NS7o3thqPqB+TBq6/h2kD6ixwIcjosrxtj4PfD0iXhxSJ48l+H/AdZLeExFPOK6XTWwNudFAxSQtAf5CGgepN5sCO0W58+F0ZWzdGhd0b2ySbicNbVPvA+pB4D2RJoqrYt6lfwHvjYhphfXDSKMyrEKaBmCDMmPr1ri6Pbbe+AqnencBsyJiXG+FJL2fkofPoHtj69a4oHtj25r0rfhFEXGNpLeTPqBukrR3ifHUepo0AOYyImK2pHcAvwZuAk51XC/q5tgach1O9W4G3t5EuaD8scu6NbZujQu6N7aGH1Ck22uPkz6gdiwxph4zgP3rbYiIecDupPHK/qfMoOjeuKC7Y2vICad6ZwBHNlHuCmCLAY6lqFtj69a4oHtj6+YPqIuALSXVm9OIiHgW+DfS1ApzHBfQ3bE15Docs5WApIOAz5LqauqOSi5pFeBbwLsjouxEbSsBJxwzMyuFb6mZmVkpnHDMzKwUTjhmZlYKJxwzMyuFE471StJ4STMk/UvSPEl/lHTWAJ3rYEnjmyg3UVLUPB6R9HNJWzV5nsm5533lmn3OuWzP8763wfZ78/aJAxVDi8dd5nXu9HkkvULSEfk9+ayk+ZLukvQ/kvrVx0nJnyR9tMH2ybk3f71t58mT6vXKCccakvQFUjv+q4ADgH8Hfklq3z8QDgbGN1n2aWBUfnwe2B64RtJaTex7agvnGWitPGeARcAWkkbWrpS0IzAsbx/oGJpVfJ07fZ6fAF8GLiG9Jz9K6t/0juh/89uDgVcDP+rHvl8HDpH0un6ee4XnoW2sN0cA346I42vWXSbp5KoCqrEkIm7Ov98saQ7wO2Bf4GfFwrmPySoR8XxE3F9inJ32DPAH4AOkjpo9PgBcC4yoIqgeZb3OkvYB3g/sGxFX1mz6RX+vbrKjgB9GxOKac61KSp4fATYBPijpfuDkiHhxeKI8rMwNwH8AR7cRwwrLVzjWm1cB/yiurP322HPbRNL+kmZJWiTpBknbFvfLt1RmSnpO0t8kfSX/MyNpMnAgsGvNrbKJLcQ6I/8cVieuu0jf/Heq3VaI7V2SrpO0QNLTkqZJ2qFm+zslTZe0UNITkr4jae3eApI0StKvJP1d0jP5Vs0hta9dP5/zFODgng/W/PPgvL5jMeTX4OLC8UbnMtvVvpZ9vc6NziNpX0lLJW1ROM8Wef3YBq/BrvnncqNz9/fqJl+ZvAO4uLDpM8CxpFEYrgA+BnwPWK/OYX5OusrxZ2sdvsKx3vwBODJfPVzey3DnmwNnkebleBY4GbhK0ut7hr2XtCfpFsiFwDHAm0nfGtcDDs+/v5aU5D6Vjzu3hViH5Z//KKw7Azglr687z4uk0cBU4DrSbZlngJ1JIzr/UdLOwNXApaRv1esBpwFD8nIjmwO/B84nfRDvDHxf0tKI+DH9f86XkEYE2IV0VfdO0qjAlwBfKymGWsPo+3VudJ6/A4+QXveJNeXHA4+RBqGs55n882uSzoyIh1qMuZ7d83H/XFi/K2mk7TPyF6nf5zHo6rkRGAoMr3Mciwg//Kj7ICWFB0gDTS4ljYR8CrBOTZnJefs7atZtDiwBDq9ZdzNwXeH4xwIvAK/JyxcD05qIayJpsMlV8+MNpGQxH9i4ENf2dfafDNxes3wT6faUGpzvd3Vi3y0ff7smX0vlWL9N+vDqWd/Uc6593vn3XwL/m3//JnBp/v1xYGInYgCmARcX1o2ufd4tvs6NzvNlUpJSTZyzSfO9NHotNgLuyOcO4E7geGBwG+/3ScBtddZ/G/hbPudkYFgvx1g1v/c/0d84VuSHL/usoYi4A3gjqUL2m6QPgi8Bt0saXFP0sYi4sWa/h0i3uN4GL97XfyvL1638hHRbd1Q/wlsPWJwffwG2BMZFxN9ryjwcEX/q7SC5kcFOwA8if2IUtq+Z4/uppFV7HsAN+dwN60wkDVFqMfVQTawTSAmyXVOA90tanXSVtdzttBJi6NHn69yH75G+pIzOy2Py8vcb7RAR/wB2APYiXe29CvgKcKOk1QAkfTTfQvxTvo07K/8+Q9Ir6xx2I1LCLvoK6crnQdL/wufzVW+9uJYAT+VjWYETjvUqIp6LiMsi4oiI2Bb4OPB64LCaYo/V2fUxYOP8+/rAK4FHC2V6luuOeNuHp0lD6Y8EXkP61nlloUzxfPUMISXSv/eyfRVSwl1c83iO9Jw26+XYk4FxpNtce+Z4vwd0YgrgXwGDSR+GawGXVRBDj2Ze54Yi4gHS1dShedWhwK0RcVcf+70QEb+NiE+Rbtd9n3Qra1Te/oOI2J70ZWcJsHNEbB8RI6KmUUCNNUh/1+J55uTjvo90xb8LcIMadw94js6+visM1+FYSyLiu5LOALapWb1hnaIbkm7BQfrWuLhOuaH5Z93Ri/uwJCL66kvTTOXxPNLtwo0bbH8qH2ciqcK46JF6OynNK/8e4NMRcX7N+o58yYuIZyRdThoB+mcR8UyxTAdiWASsVlg3pF44TR6vNxcA31Fqin8ALbbyioilkn5LSlbFD/vXA/Oi7ymXn6TBlUlOUL9Rmqp7Immqh7MlnZMTUq1X0b/39ArPVzjWkKTlEomkDYB1WfZb7YZKswz2lHkt6VvlrZC+iZJusR1UONzBpA/7m/Ly85T8zTB/UN8C/HtPq686228Gto6I2+s86iYcYHXS/9eL35hzq7ZiH6Z2nvO3SFc25zfY3m4Mc1n2iwWkq6T+6u25XpK3TyHFXPcWIYCkoQ02/RuwkPT3rPUWmqvA/wt15iiq974Abss/X10ouwGwJvDXJs630vEVjvVmpqRfAr8l3SLbnNTJciHwg5pyjwMXSTqBl1qpPUa6ndPjJFLLte+TPkyGk1oufScielpFzQLGStqf9GH3SC8f6J10HKkV2pWSJpHu148iVXhfTmrccI2kpaSK73+RbuHsB3wxIpb7cImIpyXdBpwoaT4psR5HuhW4Tk3Rfj/nSPPZT+tle7sx/AI4TNLZpNZiY4B2pqFu+FwjYpGk/wM+Dfw4Ip7q5Tg/lfQv4KekxgUbAocAY0mV9cV930JqYNCX35Neqw0i4p81638k6Y/A9aTblyNIV5YPA/cUjjGSdMV3I7a8qlst+NG9D9I//29Jt40Wkf65fwRsU1NmMqmF1wGkb3XPkf5xl2u9RapLmEn6JjuXVP+was329Ukfck+Sb2M1iGsiubVWL7FPpqaFVF/bSE1frycl06dIrd62r9m+E/AbUku4Z4C7SU3B1+0lhtcB1+Tyc0iJa5nYm33OLTzvZVqptRsD8AVSC61/kWaZ/DeWb6XW1Ovc13MF9sjr9+jjOX4s/y3m5vfSk6SEOLpB+cuADzTxfl8NeAL4SGH9+/L5/kFK2vNJiX6HOsf4BoUWjX689PAEbNaW3KFvu4gY2VdZs97kusGDgS0jYmkHjzsH2Csiilcj9cp+A3hdROzXYPtkUqKcXWfbKsBDwHERcVFbQa+gfEvNzColaWtgW9KQMCd3ONkMIXWKbbZO5WvAXyW9IercKu3DQaRbyg3rn1Z2bjRgZlX7NulW7RWk4WM6JiLmRcSgSA1Xmik/l3TLrlGrxUtJt1zrEXBYpL44VodvqZmZWSl8hWNmZqVwwjEzs1I44ZiZWSmccMzMrBROOGZmVgonHDMzK4UTjpmZleL/AW46+1xJ9j1lAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot probability distribution\n",
    "x = uncertainty_model.values\n",
    "y = uncertainty_model.probabilities\n",
    "plt.bar(x, y, width=0.2)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.grid()\n",
    "plt.xlabel('Spot Price at Maturity $S_T$ (\\$)', size=15)\n",
    "plt.ylabel('Probability ($\\%$)', size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Payoff Function\n",
    "\n",
    "The payoff function equals zero as long as the spot price at maturity $S_T$ is less than the strike price $K_1$, then increases linearly, and is bounded by $K_2$.\n",
    "The implementation uses two comparators, that flip an ancilla qubit each from $\\big|0\\rangle$ to $\\big|1\\rangle$ if $S_T \\geq K_1$ and $S_T \\leq K_2, and these ancillas are used to control the linear part of the payoff function.\n",
    "\n",
    "The linear part itself is then approximated as follows.\n",
    "We exploit the fact that $\\sin^2(y + \\pi/4) \\approx y + 1/2$ for small $|y|$.\n",
    "Thus, for a given approximation scaling factor $c_{approx} \\in [0, 1]$ and $x \\in [0, 1]$ we consider\n",
    "$$ \\sin^2( \\pi/2 * c_{approx} * ( x - 1/2 ) + \\pi/4) \\approx \\pi/2 * c_{approx} * ( x - 1/2 ) + 1/2 $$ for small $c_{approx}$.\n",
    "\n",
    "We can easily construct an operator that acts as \n",
    "$$\\big|x\\rangle \\big|0\\rangle \\mapsto \\big|x\\rangle \\left( \\cos(a*x+b) \\big|0\\rangle + \\sin(a*x+b) \\big|1\\rangle \\right),$$\n",
    "using controlled Y-rotations.\n",
    "\n",
    "Eventually, we are interested in the probability of measuring $\\big|1\\rangle$ in the last qubit, which corresponds to\n",
    "$\\sin^2(a*x+b)$.\n",
    "Together with the approximation above, this allows to approximate the values of interest.\n",
    "The smaller we choose $c_{approx}$, the better the approximation.\n",
    "However, since we are then estimating a property scaled by $c_{approx}$, the number of evaluation qubits $m$ needs to be adjusted accordingly.\n",
    "\n",
    "For more details on the approximation, we refer to:\n",
    "<a href=\"https://arxiv.org/abs/1806.06893\">Quantum Risk Analysis. Woerner, Egger. 2018.</a>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set the strike price (should be within the low and the high value of the uncertainty)\n",
    "strike_price_1 = 1.438\n",
    "strike_price_2 = 1.896\n",
    "strike_price_3 = 2.126\n",
    "strike_price_4 = strike_price_3 + strike_price_2 - strike_price_1\n",
    "\n",
    "# set the approximation scaling for the payoff function\n",
    "c_approx = 0.25\n",
    "\n",
    "# setup piecewise linear objective fcuntion\n",
    "breakpoints = [uncertainty_model.low, strike_price_1, strike_price_2, strike_price_3, strike_price_4]\n",
    "slopes = [0, 1, 0, -1, 0]\n",
    "offsets = [0, 0, strike_price_2 - strike_price_1, strike_price_2 - strike_price_1, 0]\n",
    "f_min = 0\n",
    "f_max = strike_price_2 - strike_price_1\n",
    "iron_condor_objective = PwlObjective(\n",
    "    uncertainty_model.num_target_qubits, \n",
    "    uncertainty_model.low, \n",
    "    uncertainty_model.high,\n",
    "    breakpoints,\n",
    "    slopes,\n",
    "    offsets,\n",
    "    f_min,\n",
    "    f_max,\n",
    "    c_approx\n",
    ")\n",
    "\n",
    "# construct circuit factory for payoff function\n",
    "iron_condor = UnivariateProblem(\n",
    "    uncertainty_model,\n",
    "    iron_condor_objective\n",
    ")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY8AAAE+CAYAAAB1DJw3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3dd3hUdfbH8fcBVDpiQ2zEsvYtiroWLCCK7NqwgaKoqIhSTFa36e7P7u66hVBUQBFEUex1VWyAvYCu66LYaNJEmohBWs7vj3Oj4zBJZsJkvnNnzut55knmzr3JJ8MwZ+79NlFVnHPOuUw0CB3AOedc/HjxcM45lzEvHs455zLmxcM551zGvHg455zLmBcP55xzGfPi4fKWiFwjIppwmy8iD4vIrgEzHSYi74rIdyKi0bZmIjJeRJZEOc+r5tgxSX9P1e2OnP4RP+TpIyInp9g+S0T+ESKTi49GoQM4V4uvgeOi73cBrgdeFJF9VPXbAHlGAIuALsDqaNslwAlAL2Ae8HkNx08Hzk/atijLGdPVB/gf8FjS9m7AktzHcXHixcPlu3Wq+mb0/ZsiMgd4BfgV8GCAPHsCI1V1ctK2j1X14TSO/zbh78lLqvpe6Awu//llKxc3U6OvJQAicoiIPCEiC0TkWxH5j4j0rNpZRLaILjGdl/hDxMwQkUEJ2zqJyFvR/l+KyK0i0jx67KjoMlVDYHB0uWmMiMwCLgD2q7oMVdc/rOp3iMi+SdsnichDCffHiMgUETlGRP4b/d2visg+Scc1FJE/isgnIrJaROaKyJiqnwm0B85NuHx2XvTYBpetROQMEfkg+jlfiMiNItIo4fHzop/xUxF5Pso0XUROqevz4fKbFw8XNyXR14XR13bAa9gb+AnAw8BoETkTQFWXAo8C5yX9nKOAnYE7AaI33meBxcCpwNXAWUDVm/a7wCHR9/+Mvr8eu8TzNHY56pCEfaolIo0Sb2n91RvaCfg7cCNwJrANcL+ISMI+I4BrgQeA44HLgabRY5dGmZ9OyP3vavIeC9yPPQcnAUOBK4BhKXa/F3gCe14+BcaLyA51/BtdHvPLVi7vJbzB7gLcCnwDvACgquMT9hPgZWAH4CLgvuihUcBzIrKLqs6Itp0PTFXVD6L7fwZmAyeq6vro5y3F3pAPUdU3sMtmALMSLz2JyFdAmzQvR7UH1ib9fT9R1c/SODbRFsBhqvpp9DMaYEVyD2C6iOyJFdTLVHVIwnH3A6jqhyLyLfBVGrmvAyap6rnR/Wej5+EvInKDqs5N2HeQqlYV5KnAl1jhGp7h3+fynJ95uHy3JfZmuxb4GCsg3VV1AYCItBaRISIyO2G/PsDuCT/jRawwnBsd0wI7uxidsM9BwKNVhSPyMLAO6JDFv+cj4MCk2xd1+DmzqgpH5MPoa9Wn/I7R1zF1+NnfE5GGwP5s2L50P/b+kXym9VzVN6q6BOsM4GceBcjPPFy++xroDCh2qWq+/ngq6DHAwdglpA+BFVjvp5OqdlBVFZHRQG8RuQY4A2u7uDfh57TFPiWTcNx6EVmCfcrPlgpVnZKFn7M86f6a6Gvj6OuWWOP8io38PVsBm5D03CTcT35uUuVqjCs4XjxcvltX3ZutiDTGLon0U9XhCdtTnVGPxtoxOmLtH4+p6rKExxdg7QaJP78h9ia8dGP+gAx8F33dNGl7a6wtJhNLgGYi0nIjC8hi7Gxum6TtbaKvuXpuXJ7xy1YuzjbDXsNV4y2qLkmdmLyjqn6BXVK5FrsMNTppl7eAblHBqHIK9gHr1ezGrlZV28FeVRtEZEesK3CmXoq+9qphn1rPCqLLeFOB05MeOgOoBN6oQzZXAPzMw8WWqn4tIu8A/yciK7A3sz9gl7papjhkFHbtfi7wfNJjNwDvAY+JyG3Ydfq/AROixvJ6p6pzRWQKcL2IVGCF8Urq8OleVT8WkZHAP0VkG6wjwebAaaraI9ptOtBFRLpgZyozo3aKZFcDE6JLf+OBn2KXCW9Paix3RcTPPFzcnQXMAMYCg7FG7rHV7PsU1gB+l6pWJj6gqtOArtjlmUewYnIfcFr9xK7WmcAc4B7gJqyn08d1/FmXYmdaZ2NdcsuBioTHb8Aa8B8A3sG6Om9AVZ8DegAHAE8CpVh35f51zOUKgPgytK5YiMivsAKyex26xjrnEnjxcAVPRLYDfoINbpujqscHjuRc7PllK1cM+mBjPb4DBgTO4lxB8DMP55xzGfMzD+eccxkriq66W221lZaUlNTp2G+//ZZmzZplN1A9ilPeOGWFeOWNU1aIV944ZYWNyzt16tTFqrp1ygdVteBv7du317qaOHFinY8NIU5545RVNV5545RVNV5545RVdePyAlO0mvdVv2zlnHMuY148nHPOZcyLh3POuYx58XDOOZcxLx7OOecy5sXDudqMGwclJRzZqROUlNj9fBWnrC7WimKch3N1Nm4c9OkDFRUIwOzZdh+gZ8+QyTYUp6wu9rx4OFeTq66Cioofb6uogH794OO6zpReT4YMSZ31qqu8eLis8+LhXE3mzEm9/euv4YYbcpulNtXNU1fd3+DcRvA2D+dqstNOqbe3aweVlfl1a9cus7/BuY3gxcO5mlxxxYbbmjaFG2/MfZba3HijZUvUpEl+ZnWx58XDuZrMn29ft9sOFbFP9yNH5mcbQs+elq1dO76/gHXCCfmZ1cWeFw/nqlNRASNGQLduMG8ek196CWbNyu834549YdYsJk+cCEcfDa+9BmvXhk7lCpAXD+eqc/fdsHQplJWFTlI3ZWUwbx489FDoJK4AefFwLpXKSigvh/btoUOH0GnqpmtX2H13GDSo+p5YztWRFw/nUnnuOZg+HUpLQSR0mrpp0AAuuwzeeQdefz10GldgvHg4l8qgQdC2LZxxRugkG+fcc6F1azuLci6LvHg4l2zaNDvz6N8fNt00dJqN06yZTVHyyCPW2O9clnjxcC5ZeTk0bvzDvFBx16+fXXobOjR0EldAvHg4l+irr6yXVa9esNVWodNkx447wumnwx13wDffhE7jCoQXD+cSjRgBq1dbQ3MhKSuDFStg9OjQSVyB8OLhXJXVq+GWW6BLF9h779Bpsuugg+DQQ2HwYFi/PnQaVwC8eDhX5YEHYOHC+A4KrE1pKcyYAU89FTqJKwBePJwDG0Q3aBDstRcce2zoNPWjWzebm2vQoNBJXAHw4uEcwMsvw3vvxXtQYG0aNYIBA2DyZPtbndsIXjycA+ueu+WWcM45oZPUrwsusLEfPmjQbSQvHs59/jk8/jj07WvrXxSyzTeH3r3hvvtgwYLQaVyMefFwbsgQu6Rz6aWhk+TGwIGwbh3cemvoJC7Gcl48RGRvEXlRRCpEZL6IXCciDTM4voGITBERFZHj6zOrKwJffw133gndu8N224VOkxu77WaLRA0fDqtWhU7jYiqnxUNEWgMvAAqcBFwHXA5cm8GPuRDYIfvpXFEaNQpWrrSG8mJSVgaLF8O4caGTuJjK9ZlHX6AJcIqqPq+qw7HC8RsRaVnbwVHxuRG4qn5juqKwbp1dsjr8cFu3o5gceST84hfWcO5rfbg6yHXx6ApMUNUVCdvGYwXlyDSOvx54DXixHrK5YvP44zB7duEOCqyJiJ1tTZsGzz8fOo2LoVwXjz2B6YkbVHUOUBE9Vi0R+RnQG7ii3tK54jJoEOy8M5x4YugkYfToAW3aeLddVyeiOTxlFZG1wG9VtTxp+1xgrKpeWcOxk4G3VPV3IlICzAROUNWUcy2ISB+gD0CbNm3ajx8/vk6ZV65cSfPmzet0bAhxyhsya4vp02l/ySV81q8fc087La1jCvG5bTd2LDuPHs3bY8ZQ0a5dDpKlVojPbb7YmLwdO3acqqoHpHxQVXN2A9YCpSm2zwVuquG4HsBCoGV0vwRrdD8+nd/bvn17rauJEyfW+dgQ4pQ3aNazzlJt0UL166/TPqQgn9svv1TdbDPViy+u1zy1KcjnNk9sTF5gilbzvprry1bLgFYptreOHtuAiGwC/B34G9BARDYHqhrXm4lIi/oI6grYvHk2CeKFF0LLWvtpFLZttoGzz4axY2HJktBpXIzkunhMJ6ltQ0R2BJqS1BaSoBnWNfdfWIFZBrwfPTYe8El6XGaGDYPKSpvnyVnD+apVMHJk6CQuRnJdPJ4BuiSdLXQHVgGTqzlmJdAx6XZm9NiVQM/6ieoKUkWFLfh08snWWO5g332hc2crqmvWhE7jYiLXxWM4sBp4REQ6R43a1wD/0oTuuyLymYiMAlDVdao6KfEGvBnt+oGqvpXbP8HF2tixsGxZcXbPrUlZGcyfDw89FDqJi4mcFg9VXQYcDTQEnsQGCA4Crk7atVG0j3PZU1lp3VLbt4fDDgudJr8cdxzssYd1X/ZBgy4NjXL9C1X1Q6BTLfuU1PL4LKBAF11w9WbCBPj4Y7jnnsJds6OuGjSwddsvvRReew06dAidyOU5n1XXFY9Bg2zyw9NPD50kP/XqBa1b+6BBlxYvHq44/O9/Ng1Hv36w6aah0+SnZs3g4ovh0Udh5szQaVye8+LhisPgwbbQ08UXh06S3/r1s0tYQ4eGTuLynBcPV/i++gruvtsuy2y5Zeg0+W2HHeyy3h13wIoVte/vipYXD1f4hg+H1autQdjVrqwMvvkGRo8OncTlMS8errCtXm3LrR53HOy1V+g08XDggdaVefBgWL8+dBqXp7x4uMJ2//2wcKEPCsxUaak1mj/xROgkLk958XCFS9W65+69NxxzTOg08XLyydCunXfbddXy4uEK18svw3/+Y5+ifVBgZho1goED7Tl8993QaVwe8uLhCtegQda76uyzQyeJpwsugObN7Xl0LokXD1eYPvvMrtf37WvjO1zmWrWC3r2t3Wj+/NBpXJ7x4uEK09ChdumlX7/QSeJt4EBYt856rDmXwIuHKzxffw133gk9ekDbtqHTxNuuu8KJJ9pYmVWrQqdxecSLhys8d9wBK1daQ7nbeGVltkTtPfeETuLyiBcPV1jWrbNLVkccAfvvHzpNYTjiCNhvP+u262t9uIgXD1dYHnsMZs/2QYHZJGJncR9+CM89FzqNyxNePFxhGTQIdtkFTjghdJLC0qMHbLutDxp03/Pi4QrH22/D669bD6GGvopxVm26qfVce/ZZ+Oij0GlcHvDi4QpHeTm0bGljE1z2XXwxbLaZn304wIuHKxRz58KDD9qo6BYtQqcpTFtvDeecA2PHwuLFodO4wLx4uMIwbBhUVtolK1d/Skvhu+9g5MjQSVxgXjxc/H37rb2ZdesGJSWh0xS2ffaxGYqHDYM1a0KncQF58XDxN3YsLFvmgwJzpawMFiyABx4IncQF5MXDxVtlpTXgHnCArX7n6l+XLrDnntYt2gcNFi0vHi7enn0WPvnEPg37mh250aCBrQf/7rvw6quh07hAvHi4eBs0CLbbDk47LXSS4tKrF2yxha/1UcS8eLj4+uADeOEF6N/fBrG53Gna1MZ9PPYYzJgROo0LwIuHi6/Bg22hpz59QicpTv362Uj+oUNDJ3EBePFw8bRokU0R3quXLTXrcm/77eGMM2DUKFixInQal2NePFw8DR8Oq1d799zQysrgm29s8S1XVHJePERkbxF5UUQqRGS+iFwnIjXOYici+4jIs9H+q0VkjojcISK+TFwxWr3alkXt2tW6jLpwDjgAOnSAIUNg/frQaVwO5bR4iEhr4AVAgZOA64DLgWtrObQVMBO4AugCXA10Bp4WkUb1Ftjlp/Hj4csv/awjX5SWwsyZ8PjjoZO4HMr1G29foAlwiqquAJ4XkZbANSJyc7RtA6r6OvB6wqZJIjIXeA74GfBuPed2+ULVuodWTZPhwjv5ZJsWprwcTjkldBqXI7m+bNUVmJBUJMZjBeXIDH/Wkuir99EsJpMnw/vv26ddHxSYHxo2tAkpX3kFpk4NncblSK6Lx57A9MQNqjoHqIgeq5GINBCRTUVkD+CvwDvA2/UR1OWpQYNgq62gZ8/QSVyi3r2heXMfNFhERHM4N42IrAV+q6rlSdvnAmNV9cpajn8Wa/MAmAr8SlUXVbNvH6APQJs2bdqPHz++TplXrlxJ8+bN63RsCHHKm2nWJvPmcdA55zD77LOZFWDBp0J+brNht2HD2O6xx3jzvvtYs/XWGR3rz2392Zi8HTt2nKqqB6R8UFVzdgPWAqUpts8Fbkrj+J8AvwTOxs5gpgKNazuuffv2WlcTJ06s87EhxClvxlkHDFDdZBPV+fPrJU9tCvq5zYbPP1cVUb3yyowP9ee2/mxMXmCKVvO+muvLVsuwnlPJWkeP1UhVP1XVt1T1HuwMZD/grOxGdHlp+XIbS9CjB7T1Htp5aZdd4KSTbAxORUXoNK6e5bp4TCepbUNEdgSaktQWUhtVnQ0sBXbJWjqXv+64wxZ9KisLncTVpKwMli6Fu+8OncTVs1wXj2eALiKSuMh0d2AVMDmTHxQ1mm+Jjf9whWzdOps/6cgjYb/9QqdxNTn8cNh/f+u2W1kZOo2rR7kuHsOB1cAjItI5atS+BviXJnTfFZHPRGRUwv1/iMhfRaSbiHQUkUuBCcDnWFdfV8gefRTmzPFBgXEgYv9O06fDc8+FTuPqUU6Lh6ouA44GGgJPYiPLB2EjxhM1ivapMgU4HBgF/BsYCDwMHKyq39ZzbBfaoEF2Pf2EE0Inceno3t3apbzbbkHL+dQeqvoh0KmWfUqS7o/HzzCK01tvwRtv2PTrDWucAs3li003tena//QnmDbNZgNwBcdn1XX5rbwcWraE888PncRl4uKLoXFjK/quIHnxcPnriy/gwQfhwguhRYva93f5Y6ut4JxzrNfV4sWh07h6UGvxEJFeIuKr7bjcu+UWmwhxwIDQSVxdlJbCd9/BiBGhk7h6kM6Zx2hgVwARWS8iB9VvJOewMR0jR0K3bjZjq4ufvfeGY4+1DwFr1oRO47IsneKxDNgu+l6wtTicq1933QXLlvmgwLgrK4MFC+D++0MncVmWTm+rF4C7ReRjrHCMEZFqu8eqqp+ZuI1TWWkNrQceCIceGjqN2xhdusBee1m33bPP9mn0C0g6xaM3cCmwB7A/NqL7q/oM5YrcM8/AJ5/AuHH+ZhN3InDZZdC3r633ccQRoRO5LKm1eKhqBfAPABHpDFylqu/XdzBXxAYNgu23h9NPD53EZcM558CVV9q/qxePgpFOb6v1InJgdHcSkHKpWOey4oMP4MUXoX9/2GST0GlcNjRtamcejz8On38eOo3LknQazNcAm0Xf9wIyW+XFuUyUl0OTJtCnT+gkLpsuvdRmCBg6NHQSlyXptHl8CFwjIo9hva1OE5HUK0uBquptWUvnisuiRdbOcf75sMUWodO4bNp+e5vzatQouPZaaJVqWR8XJ+kUjwHACGwCQwWuqGFfBbx4uLoZPhxWr7YGVld4ysrsw8GoUfCb34RO4zZSrZetVPV1Vf2pqm6CnXkcrKoNqrn5zHWublavhltvha5dYc89a9/fxU/79tChAwwZYmu0uFjLdG6rjthlLOey67774MsvfVBgoSsrg9mzrfHcxVpGU7Kr6mQAEfkl0AHYAlsK9lVVfSv78VxRULWG8n32gc6dQ6dx9emkk2Dnna3b7qmnhk7jNkJGxUNEmgEPAl2A9cASbCnYhiLyLHB6NC7EufRNmgTvvw+33+6DAgtdw4YwcKCdgbzzjs0i4GIp08tWNwOHAD2AxqraFmgc3T8E+Ft247miMGiQTeHds2foJC4Xeve2KfbLy0MncRsh0+JxKvB7VX1QVSsBVLVSVR8E/gD4kGCXmU8/haeegksusfEdrvC1bAkXXAAPPADz5oVO4+oo0+LRCviimse+AFpuXBxXdIYMgUaNrHi44jFwoE2AecstoZO4Osq0eLwPXCLy4wvT0f1LosedS0ujlSth9Gg480xo2zZ0HJdLO+9sjecjRkCFN5PGUabF40qssXy6iPxVRMpE5C/AR8Cx0ePO1WzcOCgp4bATTrBFn3bfPXQiF0JZGSxdCjvuyJGdOtmiX+PGhU7l0pRpV92XRGR/4M9Y+0ZbYAHwFnCKqvoYEFezceNs3qqKCr4/fb3pJnvj8Abz4jJnjvWuW7rUXguzZ/8wp5m/FvJepmceqOo0Ve2hqruqatPo61leOFxarrpqw8sUFRW23RWXq66yMT6J/LUQGxkVDxH5p4jsXV9hXBGYMyez7a5w+Wsh1jI98+gGfCAib4tIXxHxqTFdZnbaKbPtrnD5ayHWMioeqroL0BmYjq0uuEBE7o1WGHSudjfeaKOMEzVtattdcbnxRvu3T+SvhdioS5vHRFXtBWyLTde+AzBBRGaLyLUisku2Q7oCcsQR1r+/RQtUBNq1g5EjvYG0GPXsaf/27drxfcvHP//pr4WYyLh4VFHVlao6CrgaeA3YEfgj8ImIPC4i7bKU0RWSYcOsh81//8vkl16CWbP8zaKY9ewJs2bxzujRdn/x4rB5XNrqVDxEpERErhaRGcBzwEqs624L4ESgBBifrZCuQHz7rX3SPOUU65rrXKSipAS6dLER52vWhI7j0pBpb6teIvIS8BlwLjAa2FlVf6WqD6vqalV9GhgIVLdUrStWd90Fy5dDaWnoJC4flZbCwoVw//2hk7g0ZHrmMQJYCHRR1V1U9XpVnZtiv0+AG1L9ABHZW0ReFJEKEZkvIteJSI0rEIrIgSIyWkQ+i477ODrzaZxhfhdKZSUMHmxTcB96aOg0Lh916QJ77WWzLCeP/3B5J6MR5sB2qrqstp1UdQFwbfJ2EWkNvICtRngSsCvwT6yI/amGH9k92vdvwKfAz4Dro6++okwcPPMMfPKJjTD3NTtcKiK2fn3fvvDKK9a5wuWtTKcnqbVw1KIv0ASbymQF8LyItASuEZGbo22p/FVVE1vSJonId8AIEWmnqrM3Mperb+XlsN12cLrP2u9qcM45cOWVdvbhxSOvZdxgLiLdReQFEZkjIouSb7Uc3hWYkFQkxmMF5cjqDkoqHFXei75ul9Ef4HLvgw/ghRegf3/YZJPQaVw+a9rUzjwefxxmzAidxtUg0wbzs4C7sAbzHYAngKein7MCGFbLj9gTG2D4PVWdA1REj2XiEKAS+DzD41yulZfbQk9Vk945V5NLL7WBpEOGhE7iapDpmcdvsbaGftH9W1W1N7AzsBgrAjVpDSxPsX1Z9FhaRGRbrI3kblWt7WzHhbRokbVz9OoFW24ZOo2Lg+23h+7d4c47YUV1V7JdaKIZ9GoQkZXA8ao6SUTWAseo6qTosW7AIFUtqeH4tcBvVbU8aftcYKyq1roeiIhsijW67wC0r64dRkT6AH0A2rRp0378+LoNO1m5ciXNmzev07Eh5Fvednfdxc5jxvD2XXdRkTRnUb5lrU2c8sYpK2yYt8XHH9O+b18+u/RS5uZZO1ncn9tMdOzYcaqqph52oapp34D5WDddgFnAJQmPnQJ8U8vxi4CrU2z/Fisqtf1+wdpIlgB7ppu7ffv2WlcTJ06s87Eh5FXe775TbdNGtWvXlA/nVdY0xClvnLKqVpO3QwfVkhLVdetynqcmBfHcpgmYotW8r2Z62eodrHssWHvH/4nIRSJyLvB34M1ajp9OUtuGiOwINCWpLaQa5VgX35NUNZ39XUjjx8OXX/qgQFc3paU2fc3jj4dO4lLItHj8BaiabP//gLeB27CR5ouBi2s5/hmgi4i0SNjWHVgFTK7pQBH5I9AfOFtVX80wt8s1VWso33tvOOaY0GlcHJ18sk1jU15e664u99IqHiLSREROBQ4DGolIG1VdrqonAc2AzVX1l6paW9+64cBq4BER6Ry1S1wD/EsTuu9GI8lHJdw/C7gJGAvME5GDE25bZ/D3ulyZPBn+8x/79OiDAl1dNGwIAwbYgMGpU0OncUlqLR7RFOvTgAexS1N3Ax+LyLEAavNZpdUlQq1x+2igIfAkNgp9EDYzb6JG0T5Vjo2+nge8kXT7dTq/2+VYebn1rjr77NBJXJxdcAE0b26DBl1eSefM42ZsPMXhWNvEPtgAvRF1+YWq+qGqdlLVJqraVlX/rKrrk/YpUdXzEu6fp6pSzW1MXXK4evTZZ/DEEzbYq0mT0GlcnLVqZQXk/vth/vzQaVyCdIrHIcCfVPU1Vf1OVT/C2jZ2EpG29RvPxdKQIdCoEfTrV/u+ztVmwABYv96ma3d5I53i0RZIbsv4HOs2u23WE7l4W77cBnf16AFt/bOFy4Jdd4WTToIRI6CitnHILlfS7W3l8yO79IwaZYs+efdcl02lpbBkCdxzT+gkLpJu8ZiQNPnhgmj7ixlOjOgK2bp1MHSozYa6//6h07hCcsQRsN9+1hHD1/rIC+lMyb7BuhzOpfTYYzB7tveMcdknYmcf554Lzz1nC0e5oGotHqrqxcOlp7wcdt4ZTjwxdBJXiLp3h9//3l5nXjyCy3g9D+dSeucdeO01GDjQBnc5l22bbWbTtT/7LHz0Ueg0Rc+Lh8uO8nJo0QJ69w6dxBWyvn2tiAweHDpJ0fPi4TbevHnwwANw4YXQsmXoNK6Qbb21LVU7dqz1vnLBePFwG2/YMKistMFcztW3yy6DVats3IcLxouH2zgVFfaf+OSTrbHcufq27742U/Mtt8CaNaHTFC0vHm7jjB0Ly5b5oECXW6WlNtfVgw+GTlK0vHi4uqustIbL9u2hQ4fQaVwxOe442GMPHzQYkBcPV3cTJsD06b5mh8u9Bg2s7WPKFOsi7nLOi4eru/Jym/zwjDNCJ3HFqFcvaN3aVxoMxIuHq5tp02yaiP79YdNNQ6dxxahZM+jTBx59FGbODJ2m6HjxcHVTXg6NG9t/XudC6d/fLmENGxY6SdHx4uEy99VXcPfddtlgq61Cp3HFbIcd4PTT4Y474JtvQqcpKl48XOZGjIDVq63B0rnQSkthxQoYPTp0kqLixcNlZs0aG5zVpQvsvXfoNM7BQQfBoYdat/H160OnKRpePFxm7r8fFi70QYEuv5SWwowZ8OSToZMUDS8eLn2qttDTXnv5egouv3TrBjvt5N12c8iLh0vfK6/Ae+/5oECXfxo1sok5J8pEjdUAAB+TSURBVE+216ird148XPrKy2GLLeDss0MncW5DF15oYz/87CMnvHi49MyYYWuU9+0LTZuGTuPchjbfHM4/H+67DxYsCJ2m4HnxcOkZMsSWl7300tBJnKveZZfBunVw222hkxQ8Lx6uditWwJ13QvfusP32odM4V73ddoMTTrDisWpV6DQFzYuHq92oUTZ617vnujgoLYXFi+Hee0MnKWhePFzN1q+3S1YdOsABB4RO41ztjjoKfv5z61bua33UGy8ermaPPw6zZvlZh4sPEXu9TpsGL74YOk3B8uLhalZeDiUltka5c3Fx5pmwzTZ29uHqRc6Lh4jsLSIvikiFiMwXketEpGEtx2wqIn8XkVdEZJWI+LloLkydagMDBw60nlbOxcVmm1nPwKefttUuXdbltHiISGvgBUCBk4DrgMuBa2s5tClwIVABvF6fGV2C8nJo3hx69w6dxLnMXXKJLVQ2ZEjoJAUp12cefYEmwCmq+ryqDscKx29EpGV1B6nqcmALVe0CPJqbqEVu/nwYPx4uuABatQqdxrnMbbMN9OwJd90FS5eGTlNwcl08ugITVHVFwrbxWEE5sqYDVb3bRE7dcov1tBowIHQS5+qutBQqKuD220MnKTi5Lh57Aj+6AKmqc7DLUXvmOIurTkWFLfh00kmw666h0zhXdz/7GXTqBEOHwtq1odMUlEY5/n2tgeUpti+LHssaEekD9AFo06YNkyZNqtPPWblyZZ2PDSEbeds++SR7LFnCe0ceydf1+LcX43ObK3HKCvWbd8vOnfnpSy/x4fXXs6hTp43+ef7cRlQ1ZzdgLVCaYvtc4KY0f0Z/oqtY6d7at2+vdTVx4sQ6HxvCRuetrFTday/V/faz7+tR0T23ORSnrKr1nHf9etWf/ET1oIOy8poupucWmKLVvK/m+rLVMiBV62vr6DEX2nPPwUcfQVmZr9nhCkODBjZh4ttvw5tvhk5TMHJdPKaT1LYhIjtiXXG9M3Y+KC+Hbbe1SRCdKxTnnmtTtvugwazJdfF4BugiIi0StnUHVgGTc5zFJfvoI3j2WejXz/rHO1comjeHiy6Chx+G2bNDpykIuS4ew4HVwCMi0jlq1L4G+JcmdN8Vkc9EZFTigSLSVUROA34R3T8turXLXfwCN3iwjcy9+OLQSZzLvv797VLssGGhkxSEnBYPVV0GHA00BJ7EBggOAq5O2rVRtE+i24AHgQui+w9Gt471lbeoLFkCY8fCOefA1luHTuNc9u20E5x6qo35WLkydJrYy3VXXVT1Q6DG/nKqWpLONpdFI0bY4jmXXRY6iXP1p6wMHngAxoyxMxFXZz6rroM1a2xE+THHwL77hk7jXP05+GD45S/tEm1lZeg0sebFw8GDD9pcVmVloZM4V//KyuCzz+Df/w6dJNa8eBQ7Veueu8ce0KVL6DTO1b9TT4Udd/RuuxvJi0exe+01mDLFJpBr4C8HVwQaNbL2jokT4f33Q6eJLX+3KHbl5dC6tfWycq5YXHQRNG1qr39XJ148itnMmfDoozauo1mz0Gmcy53WreG88+Dee+HLL0OniSUvHsVs6FC7VNWvX+gkzuXeZZdZT8PbbgudJJa8eBSrFSvgjjvg9NNhhx1Cp3Eu93bfHX79ayse330XOk3sePEoVqNHwzffWEO5c8WqrAwWLYL77gudJHa8eBSj9ethyBA49FA46KDQaZwLp1Mn+OlPrduur3SdES8exejJJ2HGDB8U6JyInX1/8IF13XVp8+JRjMrLoV07OPnk0EmcC++ss2wyUB80mBEvHsXmvfdg8mQYMMAGSzlX7Bo3hksugaeegk8/DZ0mNrx4FJvychvTccEFte/rXLG45BJbAG3w4NBJYsOLRzFZsMB6lfTubUtyOufMttvCmWdaL8Rly0KniQUvHsXk1lth3ToYODB0EufyT2kpVFTY+CdXKy8exWLVKhg+HE44AXbbLXQa5/LPL34BRx1lMy+sWxc6Td7z4lEsxo2DxYu9e65zNSkrgy++gEceCZ0k73nxKAZVa3b8/Odw5JGh0ziXv379a9h1V59tNw1ePIrBCy/AtGn2qUokdBrn8lfDhjZh4htvwFtvhU6T17x4FIPycmjTBnr0CJ3Eufx3/vnQqpUPGqyFF49CN306PP00XHopbLZZ6DTO5b/mzeHCC+Ghh6z9w6XkxaPQDRliRaNv39BJnIuPAQOsrXDYsNBJ8pYXj0K2dCncdRf07AnbbBM6jXPx0a4dnHIKjBwJ334bOk1e8uJRyEaOtEFPvmaHc5krK4Ply+0DmNuAF49CtXatnXIffbStV+Ccy8whh8CBB9p8V5WVodPkHS8eheqhh2DePB8U6Fxdidj/n08+gWeeCZ0m73jxKESq1s1w992ha9fQaZyLr9NOg+2390GDKXjxKERvvAHvvGODnRr4P7FzdbbJJtC/vw20/eCD0Gnyir+zFKLycptyvVev0Emci78+faBJE1/rI0nOi4eI7C0iL4pIhYjMF5HrRKRhGse1EpHRIrJMRL4WkXEismUuMsfJZgsXwsMP2wu+efPQcZyLvy22gHPPhXvugUWLQqfJGzktHiLSGngBUOAk4DrgcuDaNA5/ADgKuBA4DzgQeKw+csbZ9o89Zg19/fuHjuJc4bjsMli92pY1cEDuzzz6Ak2AU1T1eVUdjhWO34hIy+oOEpFDgGOBc1X1YVV9FDgb6CAinesl6bhxUFLCkZ06QUmJ3c9n48bBTjux4/3324jyl18Onci5wrHnnvCzn8G118brPaEe38NyXTy6AhNUdUXCtvFYQalprvCuwJeq+v07oqq+DcyMHsuucePsss/s2YgqzJ5t9/P1xVKV94svELCBgfmc17m4GTcOPv4YKivj9Z5Qj+9hjbL2k9KzJ/BS4gZVnSMiFdFjT9Zw3PQU2z+KHsuuq66yN+BEFRVw3nlw001Z/3Ub7ZNPNlz5rKLC/o6ePcNkcq6QXHWVXbZKVOTvCbkuHq2B5Sm2L4seq8txu6Q6QET6AH0A2rRpw6RJk9IOeeScOaRa9ULXreOrrbdO++fkytYffpg675w5TM7g7861lStXZvTvElqc8sYpK+R/Xn9PSPXDVHN2A9YCpSm2zwVuquG454HHUmy/B3i9tt/bvn17zUi7dqo21O7Ht3btMvs5uRK3vJGJEyeGjpCROOWNU1bVGOSN2/+xLOUFpmg176u5bvNYBrRKsb119Fi2j6ubG2+Epk1/vK1pU9uej+KW17m4idv/sRzkzXXxmE5SG4WI7Ag0JXWbRrXHRaprC9k4PXvajLTt2qEiNj3zyJH5234Qt7zOxU3c/o/lIG+ui8czQBcRaZGwrTuwCphcy3HbikiHqg0icgDW3lE/M5b17AmzZjH5pZdg1qz8fZFUiVte5+Imbv/H6jlvrovHcGA18IiIdI4ata8B/qUJ3XdF5DMRGVV1X1XfAJ4DxorIKSJyMjAOeFVVX8jpX+Cccy63xUNVlwFHAw2xbrnXAoOAq5N2bRTtk6g7dnZyJzAWmAp0q8+8zjnnUst1V11U9UOgUy37lKTYthw4P7o555wLyGfVdc45lzEvHs455zImNg6ksInIV8DsOh6+FbA4i3HqW5zyxikrxCtvnLJCvPLGKStsXN52qppyCH1RFI+NISJTVPWA0DnSFae8ccoK8cobp6wQr7xxygr1l9cvWznnnMuYFw/nnHMZ8+JRu5GhA2QoTnnjlBXilTdOWSFeeeOUFeopr7d5OOecy5ifeTjnnMuYFw/nnHMZ8+LhnHMuY148nHPOZcyLh3POuYzlfFZdlx3RCoy/AgR4UFWXiMgOwBXArsAsYKSqfhAuJYjI74GnQ+dIl4g0ARqp6jcJ27YG+gN7A5XAf4BbVfXrMCmdC8+76kZERLD1QX4N7AVsgb1RLATeBMao6ifhEv5ARA4CngeaAeuApUAX4GlgPTAN2BfYFuisqq8EioqIVAKKLRd8LzBeVT8Plac2IvI08KmqXhbdPwRbrbISW0NGgPbAGqCTqk4LmHU/oImqvp6w7Tjgj/xQ6N4HrkncJ19E/+dOAPbHXiNTsA8aef2mJCItsbmiOqnqq6HzwPeZOgGbAv9W1W+jDz39sBVXZ2AfJudn7Xfm+b9TTkRP8tPYm8JC7I1he+wF/Qz25O8BXK+q14fKWUVEnsfOGrsB32ILap2MvbmdpqprRWQz4DGgsap2DJi1Evgb8FPgGCz3u1gheUBV54XKloqILAYuUNXHo/tvYs/xyVVnIyLSCngC+E5VuwTM+ibwpKreGN3vDdwBTARewgrd0cDhwKlVf1OgrK9jz+tH0f3W2Oqg7YGV0W7NsQ9qXRLP/EIQkUtreLgJ8HdgMPApgKremotcqYjIbsCLwI7RppnAsdgHzM2Bz7H3r1VAe1Wdm5VfrKpFfwPuw14EP03Yth3wLPBwdP9I7EXeOw/yLgG6JtzfBvuUeWzSfr8GFgfOWgkcFH2/BdAneqGvi26Tom1bhH5eo4wVwBEJ99ckP68Jz+23gbOuSMwGfAYMTbHfcOD9fHkdRPdHYWfMxyVsOw5YBgzKg9dBJXYWX1nNLfGx9YGzPoCdYe4W/R+7O3o/ex1oEe2zVbTPiGz9Xm8wN12BP2jCdXm107u+wMki0lZVJwM3AZcFyphMU3yffBqZV6eVqrpUVUeq6tHADsDl2Gn2cGCBiPw7aEDzPyDxTO1L7D9ksi2xQhNSZdL9dsBDKfZ7CPvkmU9OBK5T1WerNkTf3wicEizVDx4HFgEXAA1VtUHVDXs9CHBUtC15yexc6wDcqKqfqepS4E9Yu+c/NDqDU9XFQDk/fm1vFC8epgH2SSLZeuxF0iq6/xawe65C1WAKcIWINBeRBsCVwDzgEhFpCCAijYBLsTfDvKOqC1V1sKoeCuyMrWO/XeBYAH8F/iAivaPn8Ebg7yJyjIhsKiKbRe0Kf8E+8YX0CtAz4f40INXU2wdir498sjn2Ok42FWurC0pVuwHnAr8F3hGRwxIfDpOqWq2xy+1Vqv6tk9cwmoF9aMsK721lngduEJH/quoM+P6a7BDsH6Wqobw5kA89bK7CMi/DLv1UYI1lDwKfikhVg/l22KWAvKaqs7E37b/mQZZHRGQA9iltEPAx9uHhWexNQ6Jdn8DeWEK6Engt+gAxFGsov0tEtsAuBwr2uigF/hAqZIJTRaSquC0DUi0ytBV2OS44VX1ORH6GPX//FpFnsd6MQdtjUliEnXVWWQ+MwM6aE21DFrN7gzkQdXGdgJ1VzMauc+8MrAbOVNVnov1uxlbW6h4qa5Uo8/HYB4CHVXWBiGwL/A67RDEbuENV3w0YExG5Grhds9jLIxdEZEugO3AQ9klYsDe8j4CnVHVqwHjfE5FfALcBv+THxa3q+2XY5aHBYRKaqONEsjGq2jtpvxHA3qp6eG6SpSf6v3UzdkltBFZQOqrqy0GDASLyGLA0+blMsd9QYE9VPSYrv9eLh4ku95wB/BxojDU+3htdQ3Qur4nIXlgBSS50r6vq2pDZMiEiFwGfq+pLobOkEnXdHoR9QPu15kEXaBFpAzRV1Zm17PcbrOPEi1n5vV48Co+INFDVVJ/08oaINMYa9SqBz/LxDS5q89iFhDE/qjonbCrn8oM3mCcRkX1E5FQRuTC6nSoi+4TOlUxEThGRx0TkaRE5IdrWXURmAWtFZHb0KS4oETk7Gn9Qdb+RiPwV66b5X6xBf6mI5MM1eQBEpL2IPIFdH/4IeA14A5gpIvNE5DoRaRo0ZAGRSOgcqYhIk+R/axH5RfS+0D5UrrwQsn9yPt2A3lg7Qaq+3eux6T7OD50zynpGlOtVrEthBXAR1lYzChtVel+Uu0vgrB8ClyTc/2eU98/AYVg3w2uwAUxX5sFzeyzW1jUF65p9DTZQdE2U+XKsV9N/gNZ5kPd4bNzMR9Fr4YgU+/yS8GMRjiUac5Cw7WRswOg6YG30nP869HMaZWsFPBrlWgfcDjQE7kp6X3gV2Cp03jT/plOz+ToI/gflww0YEL1IbsFG424VvVAaRt93AIZFbyr98iDvO8DwhPs9o2z/TNpvNPBC4KwVwJEJ9xcBl6XY7wpgdh48t1OBu6p5jczCztYbR296twbOekz0BvZa9PqcGt3/J9El6Wi/fCge6/nxIMFu0Rvw69G//eXY2d06UgzKDJB3CDYFyQCgV/SB4WHgi6gQbo2ND5sH3BY6b5p/U1aLh7d5ACIyA3szvrmW/X4H9FXVXXKTrNocK4BTVPWF6H4rrIG0syY0NEaXs0aoarDxEyKyAOivqg9H91djZ0OTkvY7BnhCVZvkPuWPcqwCTlTV55O2t8ZG9u+jqh+JSC/gb6raNkTOKNOr2Dxc5yds64298T2P9RT8TkR+iTWcBxvMFvW2OlhV347uvwvMU9UTkvZ7GmimqkcGiJmYYxY28O726P5+WHE+X1XvStjvIuyMeecgQS3DnWnu2g4b2JiV14G3eZi2wNtp7Pc2eTCACeuGmfgCqJobaHnSfiuxwVghPYENaNw0uv8CcGaK/c7EPt2FtgjrcZfs59jzXjXOZzY/DB4NZV/gnsQNqnonNpXOwcBL0ZiPfLQv1uU12UhsosTQtuaH8V0QzWGFzROV6DNSj1fJpXOxs6Gf1nJrV90PqAsfJGjeBy4SkZe1ml5KUYPeRVgjb2izsVlTJwCo6vqoC+FHSfvtCizIcbZkf8RGQv9PRO4AngT+JiL78sNAto7AftgMq6GNBK4XkeZYoVuDjdC+CpioP4xX2QUI3fPqO2xm5R9R1anRiOgJ2GWha3KcqzqJlzm+JvWAtW/Jjw+1M7EiPDm6fzh2me1QrJ2jymGEfx18Crytqr1q2klETgPuz9Yv9eJhLsdGEH8oIo9g04dXfYpvBeyJXaPdgfwYsf0ISdMMqOpbKfbrwY9f6DmnqktF5GDszfc32ChXgEOi2xrsEsvhqvpOmJQ/UNUbo0ssfwD+r2oz1gGhNGHXtViDekj/xa67P5H8gKrOiArI08CYHOeqzgQRWRd93wr7wDA5aZ89Cf+BB2y+tcEi8lOs0J2BfRD6c/TB4n3sDKkMCD3T9ptYUatN4iDSjeZtHhER2RUbnX0cP0xtXOULrMfN3zWP16JIJiI7ActVNS+mewAQkRJ+PJDtc83PMR6bYGdujYEZ+fQcVhGRi7EpSvbTagazikgzrNdQZ7VJ/YKIZhpI9qmq3pu036Roez50Mx+IXU7dBJutYbiInIm1KVVNjDkS+H3I13DUZfgwVR1Sy35bYW12yQW7br/Xi8eGon7dVW0Fy1U19Oypzrk8EV3C3kpVvwqdJSQvHgUmOqV+F+iZD5eBJIbLukpMlvh1LiQvHgmiN43tsEspi1M8vhXwK1Udm/NwP87xqxoeboY1iv2BaDp2VX06F7lSkRgt6wrxWuI3XdG8V6er6nWBc+R8qdRsis44EpfNnYr9HcHfRMVmKz4V+/80RlWni8jPgWv54QPPMFWdkLVfGnrgSj7cgM2w6czXR7e12EjtVkn7BR9sFeWI0ypni4GTEu6/iY2IbpGwrRXWcDohD57b57FlXDfHrnUPA+Zio7c3SXi9PIP1vgr++k3jb8rq4LA6ZtgN6yVY9br8HHtTm4EV6Hewqdi/BHbIg+fsdWCvhPuto4yVUc4V/DDIsUWonFG2LtiHr4XR87oC68G4LMp3S/T/bj22nHJ2fm/of6R8uGG9apZjXXEPAAZGL+JPgZ8k7JcvxWMK1iPlfKzvduLtZ9GL+oyqbYGzxmZZ1yhHnJb43SnNW9/Qr1sCLZW6EXljs2wuNsPAg9iKh2CdKJYBo5L2uxt4M2u/N/Q/Uj7csK65/ZO2bQu8DHwFHBJty5fiIdi634uwKRN2TnisVfTC32COo0BZ3wauTrj/BdAjxX69gK/yIO+SpDeIraPn85ik/X6VB8Wj6iyztls+nIHOB85IuN8uynVK0n7nA5/kwesguXh8BZSm2C/4tDpYV+LOCfdbR/k7Je13LNYBKCu/18d5mB1JGvynqgtF5GisWr8gIj3Jj/7nqL0SRorIA8ANwAfRQi83hE2W0l+BcSLyBTCWH5Z1XYJdqqoaJJgPy7rCD0v8voqdNSUu8fuS2oDMfFni9xvgJeCOWvbrgHVDDynIUqlZlM/L5q7ix4NFq75PnuqnKTawNCu8eJj5wE+wM43vqfXd7iEig7HTwqAN5clUdTnQX0RGYn3PPwX+Rh6tsazxWtYV4rXE79tYu9y/a9opWjsltCBLpW6kuCyb+xrwfyLyaZTlH9hs1r+PZs34Jpr/7ndYscsK723F9xOL7aKqR9Wwzx+xT82qASeYq4mI9MCWytwBmwAt+BKZVSQmy7pCrJb4/TPQR1WTB7Um73cEcK2qdsxNspQZgiyVWlcSo2VzRWQ3bCqdqtfBLOxs/iFsxP5soAT7MNRRVf+Tld/rxeP7bm7dgb9oDcvOishZ2LXv86vbJ7TokkozYKWqrg+dxzkIt1RqfZM8WTY3Gt91GNZD8EVVXRUNdr6QHz7w3Kuqc7P2O714OOecy1Q+zF7p6omI3C4io0LnSEecskL88jqXbd5gngERuR1ooKoXhM6Spo7E5wNCnLJCjPKKyAvYVYajQ2epTZyyQrzyZjurF4/MxOYNA0BVdwudIV1xygqxyyvE53Ubp6wQr7xZzeptHgUs6qK5jaqGXqymVnHKCvHL61y2xaVi5gURaRytkREXv8ZWRIuDOGWFGOUVkU3i8rqNU1aIV95sZ/XikZnYvGG44iAi/UTkcxFZJSLvi8g5KXbbnzx43cYpK8Qrb4is3uYRQyKSbp/yVCNicypOWSFeeaNBoUOxJXLfw5YiHSMiJwFnq2rWpqLYWHHKCvHKGyqrFw/i9YYROQKb5uPDWvbLh2kp4pQV4pX3CuAfqvr9vFXRfGzjgIkicryqLgmW7sfilBXilTdIVm8wB0RkHem9YWwP/DL09CQi8j4wXVW717LfacD9IfPGKWuUIzZ5ReQb4ARVnZS0vQRbb6QhNv/W1sDrnjV9ccobKqu3eZhpwP9U9fSabsC/QgeNvAkcnMZ+iRMPhhKnrBCvvF9jE/P9iKrOwi5dLAbeAA7MbayU4pQV4pU3SFY/8+D7yc2OU9V2tex3KramddCiKyK7Avuo6hO17NcE606aPO11zsQpa5QjNnlF5HHgG1U9u5rHm2CT43Ul8ISeccoa5YlN3lBZvXgQrzcM56qIyOlAGXB8dRN6ikhD4DZsQs+dc5kvKUdsskZZYpM3VFYvHs455zLmbR7OOecy5sXDOedcxrx4uKIiIueJyFQR+UZElonIeyJSL73oRGR3EblGRDZPY99rREQTbvNF5OGoPa62Y8+LjmmeneTO1c6LhysaYksJ3wFMAE4BegGPAyfW06/cHbgaqLV4RL4GDoluVwC/AF4UkWa1HPfv6JiKOuZ0LmM+wtwVk/7ACFW9MmHbkyJybahASdap6pvR92+KyBzgFeBXwIPJO0c9aBqq6lfAV7mL6ZyfebjisjmwMHmjJnQ5FJGS6BLQWSJyd3R5a5GIXJ18nIh0EpG3ROQ7EflSRG6tunQkIkcBT0a7zox+5qwM806NvpZEP3OMiEwRkZNFZBrwHfDLVJetRKSJiNwsIrNFZLWIzBSRvyTlv1BEpkWPzxaR3+FcmvzMwxWTd4EB0Sf6p2qZ7+fvwFPAadh8V1eLyGJVvQVARPYBngWeB04FdgT+CuyCTQXxLtGcQ9glsgXA6gzzlkRfFyZtuxm4Lto+E/hRu4iICHY57hDgeqwIbQ8cnrDPb4Gbop81CWgPXC8iFao6LMOcrhipqt/8VhQ34GfADGxqkUpsWprrgJYJ+5REjz+XdOztwDxsGWKA8cCn2GWjqn3OiI49JLp/fHS/JI1s12DTSDSKbrsDE4EVQNtonzHRz/tF0rHnRdubR/e7RPdPrOZ3tQRWAlcnba8qSA1ry+s3v/llK1c0VPW/wF5YA/mt2NxUfwampOip9GjS/UeA7YAdovsHAY+q6vqEfR4G1gEd6hhxS2BtdPsYO4vprqoLEvaZp6r/qeXndAKWavUzJhwCNAMeFJFGVTfgJaANP/yNzlXLL1u5oqKqq7G2iCcBROQCrAfWBcDghF0XJR1adb8tMCf6+mXSz14vIkuALeoY72ugM3bWsBCYr6rJU0B8ucFRG9oSu0xWnapJ9KZV8/iOgE/B42rkxcMVNVUdJSI3A3smPbRNNfcXJHz90T5R76ctgZTzC6VhnapOqWWfdOYTWoIVt+pU5Tue1MXo4zR+hytyftnKFQ0RSS4IiMjWQCs2fBPtlnS/qtF7bnT/LaBbVDAS92kEvBrdXxN9zfXCUS8CW4jI8dU8/gawCthOVaekuH2Tu6gurvzMwxWTD6Lpq5/DLkO1w3pEVQB3Je27TzRV/8NYb6sLgMtUtTJ6/AZsyc/HROQ2rJ3gb8AEVX0j2qfqE/zFIjIeqFDVD+rnT/uR57GBkPeKyHVYz6+2wBGqerGqLheRa4DBItIOeBn7ILk70FFVkwuncxvw4uGKyXXAScAQrF1iIfA61ig9M2nf32GXdR7GxlNcD3zfhVVVp4lIV6y76yNYr6j7ouOq9pktIlcAA4EB2FlLSX38YYlUVUWkW5S5FFtBbj5wb8I+N4vIfGwq78uxv/ET4P76zucKg0/J7lyCaOnOmdiynk+FTeNc/vI2D+eccxnz4uGccy5jftnKOedcxvzMwznnXMa8eDjnnMuYFw/nnHMZ8+LhnHMuY148nHPOZez/AV2UWkzpP7TSAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot exact payoff function (evaluated on the grid of the uncertainty model)\n",
    "x = uncertainty_model.values\n",
    "def payoff(x):\n",
    "    if x <= strike_price_1:\n",
    "        return 0\n",
    "    elif x < strike_price_2:\n",
    "        return x - strike_price_1\n",
    "    elif x < strike_price_3:\n",
    "        return strike_price_2 - strike_price_1\n",
    "    elif x < strike_price_4:\n",
    "        return strike_price_2 - strike_price_1 + strike_price_3 - x\n",
    "    else:\n",
    "        return 0\n",
    "y = [payoff(x_) for x_ in x]\n",
    "plt.plot(x, y, 'ro-')\n",
    "plt.grid()\n",
    "plt.title('Payoff Function', size=15)\n",
    "plt.xlabel('Spot Price', size=15)\n",
    "plt.ylabel('Payoff', size=15)\n",
    "plt.xticks(x, size=15, rotation=90)\n",
    "plt.yticks(size=15)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "exact expected value:\t0.3569\n"
     ]
    }
   ],
   "source": [
    "# evaluate exact expected value (normalized to the [0, 1] interval)\n",
    "exact_value = np.dot(uncertainty_model.probabilities, y)\n",
    "print('exact expected value:\\t%.4f' % exact_value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Evaluate Expected Payoff"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# set number of evaluation qubits (=log(samples))\n",
    "m = 5\n",
    "\n",
    "# construct amplitude estimation \n",
    "ae = AmplitudeEstimation(m, iron_condor)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "# result = ae.run(quantum_instance=BasicAer.get_backend('qasm_simulator'), shots=100)\n",
    "result = ae.run(quantum_instance=BasicAer.get_backend('statevector_simulator'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Exact value:    \t0.3569\n",
      "Estimated value:\t0.3428\n",
      "Probability:    \t0.9697\n"
     ]
    }
   ],
   "source": [
    "print('Exact value:    \\t%.4f' % exact_value)\n",
    "print('Estimated value:\\t%.4f' % result['estimation'])\n",
    "print('Probability:    \\t%.4f' % result['max_probability'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAAEPCAYAAAB/WNKuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAdzElEQVR4nO3dfbQcVZ3u8e/DixAGiJF3kCGAaAauzNzhgGRdhETew1yjAQwXvLPioFGvCjMLGV5ECOiwBK+ALnQFlk643NFERhjuIIQQICe8KgYB0SRg0IACgswcEkMgvOR3/9gVKCp9uqv7dNVJnzyftXp1965du3btVPp3qnbtXYoIzMzMum2T4a6AmZmNTA4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxa0LSDEn92ed+STPaXH+CpCiWNUjemyQ92mT5lZJelLRFyW2/R1JIOqadOpt1iwOM2YZjNvBfJO1bXCBpU+AE4IaIWFN7zcw64ABjtuH4f8Bq4H80WDYR2IkUhMx6ggOMWYckjZf075KelfSSpIclndJpeRHxEnATMLXB4pOA54E7s23vJmmWpN9KelnS45IulLR5k/pull0y+0wh/auS/lBI20PSDyUNSFotaa6kfTrdN9s4bTbcFTDbkEXEjNznCYXFewD3AjOBV4D/BsyStDYiZmfr9AMqltXEbGCqpAMi4kGALGhMAb4fEW9k+XYAXgD+HngRGAdcAGwPfK7N3XwbSdtn+/UcMD3bt3OB+ZLe50t0VpYDjFmHImLOus+SBNwFvBv4FJ1fyppLChgnAQ9maUcDY/JlRsTDwMO57d8LvAzMlHR6RLze4fYBzgC2AA6PiBez8u8DlgPTgKuGULZtRHyJzKxDksZI+pakJ4HXstd04L2dlhkRrwI3AB/LghakS2ZPAvfntr2JpDMkLZH0crbt/wOMIgW5oTgCmAesyi6rbQasAH4O9A2xbNuIOMCYde4a0o//14GjgAOBfwa2HGK5s4E/B8ZL2hKYDMyJt099fgZwCfCvwIeBg4DTsmVD3f72wCm8FTTXvQ4Fdh9i2bYR8SUysw5kP/x/A3wuImbm0rvxR9sCUv/HScAuwDasf8ntRFLQOT+37f1blPsG8DrwjkL6mML3/wQeAi5uUMbKFtswe5MDjFlntiBdAXizw1vSNqSziSE9ZCki3pB0HSmI7AYsiYhHCtlG5bedaXoHW0SEpKeBv8jVeVPg8ELWO0hnTY+6Q9+GwgHGrAMRsULSz4DzJa0E1gJnk/oqtu3CJmYDXwA+Sro7rGg+8FlJi4DfAH8LjC1R7r8B0yU9QurX+RSwVSHP/wZOBu6UdCXwDLAzcBjQHxHXtb03tlFyH4xZ504m/bhfC3wTuD77PGQRcT/pri3R+I60C4DrSJexZgMvAf9QoujzSTcRXAzMAhZRqHNEPA8cDCwDrgBuI/X3bAMMOpWNWZHqfmSypPcAZwLjgf2AuxuML2i03mjSwf4RUmD8MXBaRPxHId9k4KvAPqT//BdGxA+7uQ9mZtbacJzB7AdMAh4DHm9jveuACcAnSffiHwjcmM8g6RDSX5ELgGOBm4HZko4aaqXNzKw9w3EGs0lErM0+/wjYvtUZjKTxwH3AYRFxV5Z2EPBT4MiIuD1LmwdsHhEfyq17C7BtRBxSxf6YmVljtZ/BrAsubToWeG5dcMnKeQD4bbaMbArziaQznbw5pPEEozursZmZdaJXOvnHAUsbpC/JlgHsDWzeIN8S0n52PLrazMza1yu3KY8hzc9UNADslctDg3wDheVvI2k6aXoPRo0adcDuu3dnoPLatWvZZJNeid/Dx+1UjtupHLdTOd1sp8cff/yFiNih0bJeCTCViYirgasB+vr6YtGiRV0pt7+/nwkTJnSlrJHM7VSO26kct1M53WynbC6+hnol1A8AjfpQxvDWGcq692K+MYXlZmZWg14JMEt5q68lL9838wRpQr5ivnGkUdbt3BJtZmZD1CsBZi6wczbOBQBJfaT+l7kA2ZxJC0jzN+VNBe6PiBU11dXMzBiGPhhJW5EGWkKayG9bSSdk32+JiNWSlgELI+JUSNNmSLoNuFbSF0lnJJcA96wbA5P5CtAv6QrSIMxJ2euYynfMzMzeZjg6+XckPcMib933PUnzL20GbFrIMxW4nPS8jTenislniIh7smD1VeCzpHEyJ0fEbV2sv5mZlVB7gImI5WTPKG+SZ2yDtBeBT2SvZuveSGEKGTNrbuzZN7fMs/xrx9VQExtJeqUPxszMeowDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSUcYMzMrBIOMGZmVgkHGDMzq4QDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVWi9gAjaV9Jd0haLekZSRdJ2rTFOjMkxSCvc3L5rhkkz7jq98zMzPI2q3NjksYAtwOLgcnA3sA3SIHuvCarfhe4tZD2EeAsYG4hfSnwiULa8s5qbGZmnao1wACfAUYBUyJiJTBf0rbADEmXZmnriYjfA7/Pp0n6MrA0Ih4uZH8pIn5SQd3NzKwNdV8iOxaYVwgkc0hB57CyhUjaDjgSmN3d6pmZWbfUHWDGkS5hvSkingJWZ8vKOh7YnMYBZl9JKyWtkXSPpNKBy8zMukcRUd/GpNeAMyPiikL674FrI+LckuXcCYyOiAMK6acDr5L6eHYAzgAOAA6JiAcGKWs6MB1gp512OmDOnDnt7dQgVq1axdZbb92VskYyt1M5VbfTo0+vaJnn/buNrmz73eLjqZxuttPEiRMfjIi+Rsvq7oMZMkm7kC6nnVVcFhHfLOS9BfgVcC7ppoD1RMTVwNUAfX19MWHChK7Us7+/n26VNZK5ncqpup2mnX1zyzzLT6lu+93i46mcutqp7ktkA0CjP4PGZMvK+Bgg4IetMkbEauAW4K/LVtDMzLqj7gCzlEJfi6Tdga0o9M00cRJwT0T8rmT+yF5mZlajugPMXOBoSdvk0qYCLwMLW60saSxwMCXvHpM0CjgOeLDdipqZ2dDUHWBmAmuAGyQdkXWwzwAuy9+6LGmZpO81WP8k4HXgX4sLJI2WdLekT0s6XNJUYAGwK3BxBftiZmZN1NrJHxEDkg4HrgRuAl4ELicFmWK9Gk0fcxJwR0S80GDZGuCPpBkBdgReAe4HDouIRV3ZATMzK632u8giYjHwoRZ5xg6S/ldN1nkFmDKkypmZWdd4NmUzM6uEA4yZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSUcYMzMrBIOMGZmVgkHGDMzq0RbAUZSo+lbzMzM1tPuGczTki6V9BeV1MbMzEaMdgPMTOAE4JeSfippuqRtK6iXmZn1uLYCTETMiIi9gCOBx4DLgGclfV/SEVVU0MzMelNHnfwRcWdE/C2wM/AF4H3APEnLJc2QtGs3K2lmZr1nqHeR9QGHkh6DPADcDXwSWCbp40Ms28zMeljbAUbSHpIukPQEcAewC/B3wK4R8T+BPYCrgK93taZmZtZT2nrgmKQFwAeBp4FZwKyIeDKfJyLekPQD4PSu1dLMzHpOu0+0fB6YBMyPiGiS72Fgz45rZWZmPa/dS2TfBu5rFFwkbS3pUICIeK14ZmNmZhuXdgPMAmDfQZa9L1tuZmbWdoBRk2VbA6uHUBczMxtBWvbBZJe9JuSSPinpmEK2LYHjgEe7VzUzM+tlZTr5P0AaTAkQwInA64U8rwJLgTO7VzUzM+tlLQNMRHydbEyLpN8CH42Ih6uumJmZ9ba2blOOCN96bGZmpZTpg5kE3BMRK7PPTUXELV2pmZmZ9bQyZzA/Bg4GHsg+B4PfTRaAH0pmZmalAsyewLO5z2ZmZi2V6eR/stFnMzOzZsr0wWzVToER4cGWZmZW6hLZKlLfSlnugzEzs1IB5u9oL8CYmZmV6oO5poZ6mJnZCDPURyabmZk1VKaT/wFgWkQslvQzWlwui4iDulU5MzPrXWX6YH4FvJz77P4YMzNrqUwfzCdyn6dVWhszMxsxOu6DUbKDpGYPITMzs41U2wFG0iRJ9wGvAH8AXpF0n6Tjul47MzPrWW0FGEmfBm4iDb48nfTwsdOz7/+eLTczM2vveTDAucBVEfG/CukzJc0EvgRc1ZWamZlZT2v3Etl2wL8Nsux64F2tCpC0r6Q7JK2W9IykiyQ1nV5G0lhJ0eA1p0HeyZIelfSKpMWSppbaMzMz66p2z2AWAIcB8xssOwy4q9nKksYAtwOLgcnA3sA3SIHuvBLb/yJwb+77C4XyDyEFuu8ApwGTgNmSBiLithLlm5lZl5QZaLlv7uu3gO9K2g64EXge2BH4KHAs8MkWxX0GGAVMiYiVwHxJ2wIzJF2apTXzWET8pMnyLwN3RcRp2fcFkvYDzgccYMzMalTmDOaXvH1wpYBPZ6/i0y1vpflsyscC8wqBZA5wCekM6KYS9WlI0hbARNKZS94cYJak0RGxotPyzcysPWUCzMQubm8ccGc+ISKekrQ6W9YqwMyS9C7SmdNs4EsRsW6Wgb2BzYGlhXWWkC7BvRf42dCqb2ZmZZUZyb+wi9sbA7zYIH0gWzaYNcC3SZe5VgITgLNIQWVyrmwalD9QWP42kqYD0wF22mkn+vv7m9W/tFWrVnWtrJHM7VRO1e10xvtfb5mnF/6dfDyVU1c7tdvJ/yZJmwBbFtOreKJlRDwLfD6X1C/pOeA7kv4yIh4ZQtlXA1cD9PX1xYQJE4ZU1zcr2N9Pt8oaydxO5VTdTtPOvrllnuWnVLf9bvHxVE5d7dTuQEtJOkvSMuA14E8NXs0MAKMbpI/hrTONsn6UvR+QK5sG5Y8pLDczsxq0Ow7mNOBs4Hukzv1/Ai4CHgeWk11qamIpqa/lTZJ2B7Zi/b6TVqLw/gQp6I0r5BsHrM3qaGZmNWk3wHwKuAC4NPt+Y0RcCOxHChD7tFh/LnC0pG1yaVNJjwNot6/nhOz9QYCIWEMap3NiId9U4H7fQWZmVq92+2D2BB6OiDckvQa8EyAi1kr6DvBd0hnOYGaSzoJukHQJsBcwA7gsf+tydgluYUScmn2fAWxDGmS5EjgUOBO4ISJ+kSv/K6T+mStI43QmZa9j2txPMzMbonbPYP4D2Dr7/BTwX3PLxpAGUQ4qIgaAw0ljZW4CLgQuJ50V5W3G28fTLCWNk5kF3AKcDHw9e8+Xfw/pzOYIYB7wYeBkj+I3M6tfu2cw9wIHkn7kf0Aagf8u4FXgc8AdrQqIiMXAh1rkGVv4Poc0YLKliLiRdPZiZmbDqN0AMwPYLft8MekS2TTSmct84AvdqpiZmfW2tgJMRDwGPJZ9XkN6FszpFdTLzMx63FAGWr4b2AV4JiKe7l6VzMxsJOjkkcmflfQ74Engp8BTkn4vqfgQMjMz24i1O5L/fOBK0niW44C+7H0u8K1suZmZWduXyD4HXBwRXy6k35rNDfY50sh+MzPbyLV7iWwUgz+1ciENJr80M7ONU7sB5kZgyiDLjgd+PLTqmJnZSFHmkcmTcl/nApdKGsv6j0zeD/jH7lfRzMx6UZk+mB+z/qORdwOObpD3X0hPmjQzs41cmQCzZ+W1MDOzEafMI5OfrKMiZmY2srQ9kl/SZqQO/UOAdwH/CdxNmjq/9YO9zcxso9BWgJG0I3AbsD/pCZbPAeNJ418ekXRURPyx25U0M7Pe0+5typcB2wEHR8ReETE+IvYCPpClX9btCpqZWW9qN8BMAs6KiAfyiRHxM+Ac0rQxZmZmbQeYLYA/DbLsT8A7hlYdMzMbKdoNMD8BzpL0Z/nE7PtZ2XIzM7O27yI7A1gA/E7SbaRO/h1Jgy4FTOhq7czMrGe1dQYTEQ8D+wBXAzsAR5ICzExgn4h4pOs1NDOznlT6DEbS5sBBwG8j4uzqqmRmZiNBO2cwbwB3AuMqqouZmY0gpQNMRKwFfg3sXF11zMxspGj3LrIvAedLen8VlTEzs5Gj3bvIziON2H9Y0tOku8ginyEiDupS3czMrIe1G2B+mb3MzMyaKhVgJI0iTRPzS+APwO0R8VyVFTMzs95W5pHJewG3A2NzySslfSwibquqYmZm1tvKdPJfCqwFPghsBewHPARcVWG9zMysx5UJMOOB8yLi3oh4JSKWAJ8G/lzSLtVWz8zMelWZALML8JtC2hOkucc8JsbMzBoqOw4mWmcxMzN7S9nblOdJer1B+h3F9IjYcejVMjOzXlcmwFxYeS3MzGzEaRlgIsIBxszM2tbuXGRmZmalOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVWi9gAjaV9Jd0haLekZSRdJ2rTFOgdKmiVpWbbeY5IukLRlId8MSdHgdUy1e2VmZkXtPnBsSCSNIU39vxiYDOwNfIMU6M5rsurULO8lwK+B/YGvZO/HF/KuAIoBZclQ625mZu2pNcAAnwFGAVMiYiUwX9K2wAxJl2ZpjXwtIl7Ife+X9ApwlaQ9IuLJ3LLXI+In1VTfzMzKqvsS2bHAvEIgmUMKOocNtlIhuKzzUPa+a/eqZ2Zm3VJ3gBkHLM0nRMRTwOpsWTvGkx6E9kQh/Z2SXpD0mqSHJE3puLZmZtYxRdQ3E7+k14AzI+KKQvrvgWsj4tyS5ewM/AK4JSKm5dI/DuxIOrvZhvRgtEnA8RFxwyBlTQemA+y0004HzJkzp93damjVqlVsvfXWXSlrJHM7lVN1Oz369IqWed6/2+jKtt8tPp7K6WY7TZw48cGI6Gu0rOcCjKR3kG4UeDdwQEQMNMkr4D5gVET8Vauy+/r6YtGiRa2yldLf38+ECRO6UtZI5nYqp+p2Gnv2zS3zLP/acZVtv1t8PJXTzXaSNGiAqfsS2QDQ6M+gMdmyprKAcS2wHzCpWXABiBQ9bwD2b3UrtJmZdVfdd5EtpdDXIml3YCsKfTODuIJ0e/OREVEmP6SncfqJnGZmNav7DGYucLSkbXJpU4GXgYXNVpR0DvB54OMRcU+ZjWVnPMcDj0TEG51V2czMOlH3GcxM4DTgBkmXAHsBM4DL8rcuS1oGLIyIU7PvJwMXA9cAT0s6OFfmExHxxyzfQuB60tnQnwGfAj4AfKTa3TIzs6JaA0xEDEg6HLgSuAl4EbicFGSK9cr3mRyVvU/LXnmfIAUegGXA3wO7kG5h/jlwXETM7Ub9zcysvLrPYIiIxcCHWuQZW/g+jfUDS6P1Th1C1czMrIs8m7KZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSUcYMzMrBIOMGZmVgkHGDMzq4QDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwq4QBjZmaVcIAxM7NKOMCYmVklHGDMzKwSDjBmZlYJBxgzM6uEA4yZmVXCAcbMzCrhAGNmZpVwgDEzs0o4wJiZWSUcYMzMrBIOMGZmVgkHGDMzq4QDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSjjAmJlZJRxgzMysEg4wZmZWCQcYMzOrhAOMmZlVwgHGzMwqUXuAkbSvpDskrZb0jKSLJG1aYr3RkmZJGpC0QtL3JW3XIN9kSY9KekXSYklTq9kTMzNrptYAI2kMcDsQwGTgIuAM4MISq18HTAA+CUwDDgRuLJR/CHA9sAA4FrgZmC3pqK7sgJmZlbZZzdv7DDAKmBIRK4H5krYFZki6NEtbj6TxwFHAYRFxV5b2NPBTSUdExO1Z1i8Dd0XEadn3BZL2A84Hbqtut8ysHWPPvrllnuVfO66GmliV6g4wxwLzCoFkDnAJcBhwU5P1nlsXXAAi4gFJv82W3S5pC2AicFph3TnALEmjI2JFl/bDrFLFH+Az3v8603Jp/vG1XlB3gBkH3JlPiIinJK3Olg0WYMYBSxukL8mWAewNbN4g3xLSpcD3Aj/rrNq2MWv117Z/7Ddc/rcbXnUHmDHAiw3SB7Jlnay3Vy4PDfINFJa/jaTpwPTs6ypJjzWpRzu2B17oUlkjWc+3ky6pfhunFdqpjm0W1b3NDrfX1vE0HO24gejm/7s9BltQd4DZ4ETE1cDV3S5X0qKI6Ot2uSON26kct1M5bqdy6mqnum9THgBGN0gfw1tnGp2ut+69mG9MYbmZmdWg7gCzlLf6TACQtDuwFY37WAZdL5Pvm3kCeK1BvnHAWuDxDuprZmYdqjvAzAWOlrRNLm0q8DKwsMV6O2fjXACQ1Efqf5kLEBFrSONfTiysOxW4fxjuIOv6ZbcRyu1UjtupHLdTObW0kyKiju2kjaWBlouBX5JuTd4LuAy4IiLOy+VbBiyMiFNzafOAfYAvks5ILgGej4gP5vIcAvQDV5IGYU7K8h8TER4HY2ZWo1rPYCJiADgc2JR0S/KFwOXABYWsm2V58qaSznL+GbgWeBD4aKH8e4ATgCOAecCHgZMdXMzM6lfrGYyZmW08PJtyC56cs7VO2kjSgVn7LMvWe0zSBZK2LOSbISkavI6pdq+6r8N2GjvI/s9pkLfnjyXouJ0GO05C0jm5fNcMkqfRTUQbNEnvkXSVpF9IekNSf8n1avtt2ujHwTSTm5xzMWlyzr2Bb5AC83lNVoU0Oed7SZNzruszuhEo9hldD3yHNMXNJNLknAO9cllvCG00Nct7CfBrYH/gK9n78YW8K4BiQFky1LrXaYjHEqS+xHtz3982SG4kHEswpHb6LnBrIe0jwFlkNwLlLAU+UUhb3lmNh9V+pH/nn5BmMSmrvt+miPBrkBdwDmn8zLa5tH8EVufTGqw3njRj9KG5tIOytCNyafOAOwvr3gLcM9z7XkMbbd8gbXrWRnvk0mYALwz3fg5jO43N2uRvWpTf88fSUNppkLJuBpYU0q4BFg33fnaprTbJff4R0F9inVp/m3yJrLnBJuccRZqcs9l6603OCaybnJPc5JzXFdadA4yX1Ghg6YaoozaKiEbTVDyUve/aveptMDo9lloaQccSdKmdsks+RwKzu1u9DUdErO1gtVp/mxxgmltvks2IeIr011Sza7bdmpyzF3TaRo2MJ52yP1FIf6ekFyS9JukhSVM6ru3wGWo7zcqusz8r6TJJo3LLRsqxBN07no4ntUmjALOvpJWS1ki6R9KQAnyPqfW3yQGmuSom5xyTy0ODfE0n59wAddpGbyNpZ9I19v8bEc/nFi0jXSI5kfSj8QxwfQ8GmU7baQ3wbeBU0i3+VwGfJf01mS+bBuX32rEEXTqegJOAn0fErwvpD5EecvjfgVNIwyHmSzqog7r2olp/m9zJb8NO0jtIp+OrgH/IL4uIfynkvQm4j/QQuRvqquNwiYhngc/nkvolPQd8R9JfRsQjw1S1DZakXUiX084qLouIbxby3gL8CjiXdFOAdZHPYJrz5JytddpGAEgSaeDsfsCkSINxBxWpt/EGYP8yt4tvQIbUTgU/yt4PyJVNg/J77ViC7rTTxwABP2yVMSJWkzqv/7psBXtcrb9NDjDNeXLO1jpto3WuIN2OOjkiyuSHdMdLr40QHmo75UXhfaQcS9CddjqJdLfT70rm78XjqVO1/jY5wDS3MU3O2alO24hsANzngY9HmuanpeyM53jgkYh4o7MqD4uO26mBE7L3B2FEHUswxHaSNBY4mJJ3j2U3SxxH1pYbgXp/m4b7Xu4N+UU6JXwWmE+a32w6qZ/gq4V8y4DvFdLmAb8BppCu7T4G3F3IcwjwOumv+AnApaS/EI4a7n2vuo2Ak0l/Nc4i/SDkXzvk8i0kDfQ6ijT33C1ZG314uPe9pnaaQRpoOCVb7yLSj+31I+1YGko75dLPJv313Wic1WjgbuDTpBsmppIGKa4B+oZ73ztoq61If2ycANxP6kta932rwdqpzt+mYW+kDf0F7Avcmf2nfpY02nzTQp7lwDWFtHdmP54vAiuBHwxy0H+ENLv0GtIp6knDvc91tBFpwFsM8pqWy/e97D/Dy8BL2Q/EscO9zzW200nAItJsBq9mPxgXAVuMxGOp03bKpT8M3DpIuVuS+u9+l7XRCtLo/4OHe587bKexTf4PjR2sner8bfJkl2ZmVgn3wZiZWSUcYMzMrBIOMGZmVgkHGDMzq4QDjJmZVcIBxszMKuEAY2ZmlXCAMTOzSvx/VaxREbm57+YAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZgAAAEPCAYAAAB/WNKuAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjEsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8QZhcZAAAgAElEQVR4nO3deZxcVZ338c+XRVkCISwJ4CARFDMsIw5he4aRREC2YBBE1tEIGBxZfSETQMQgGg3INmF8gMeR4CjEBWWGAAlbEsSNRcOELZhoWE3YAiEm7L/nj3ObFDfV1beq61Z1db7v16te1XXuubd+p7q7fnXPueeUIgIzM7NmW63dAZiZWf/kBGNmZqVwgjEzs1I4wZiZWSmcYMzMrBROMGZmVgonGGsqSeMlRTe3YwoeY5vsOBvkysdkxxlQTvTF4ujlMX8uaWaBemtIOk3SA5KWS1os6WZJezT4vH3lNR2T+5t4VtJ0Sf9YYN8R2T7btyJW6z0nGCvDy8DuVW7TCu6/DfB1IP/GflN2nGXNCbPhOEolaXXgBmAC8D/AAcAY4C1gpqSjGjhsX3lNu3w8e94TgE2AGZI272GfP2T7zC85NmuSNdodgPVLb0bE75p90Ih4Dniu2cftg04GDgT2j4jKpPzfkqYAV0maFRFP9/aJ2via3hsRSwEk3Qc8DhwNXJivKEnAeyNiCdD0vysrj89grC0knSVpnqRXJS2SNE3SppJGADdm1f6SdYksyPZ5V3eOpKHZ4yMkXS1piaSnurriJP2bpGckPSdpoqTVKp5/mKQpkp6UtEzSQ1mX1GrZ9m7jyLa/P9v/xWz/6ZI+nGvjFlm31nJJCyQdX/DlORWYkUsuXb4KrAUcV/E8CyR9V9LXJC2UtFTSjyUN7Kkt1brIJG0s6RpJL2RtmylpeK5tXc/55ew1X5y9HnWf7UXEk6QkNzQ79nhJz0vaQ9K9wKvAYdW6yCStnv0tPSbptSyWyblYR0u6L/tbWyjpAklr1hun1c9nMFYKSSv9bUXEm9m2zwJnA+OAh4CNSF0m65K6Qb4CfBc4BPgr8FoPTzcR+DFwKHAscI2kjwJbZo93Ar4J/BGYku3zPmButt8rwI7AecDawLdrxSFpQ+Bu4AXgi6TupTOB2yVtExHLs0/d/w1sTEoGr2bH3xD4U43XbQvSG+0l1bZHxHxJc4CP5TYdCcwDvgBsBlwAfB84rFZbunED8MFsn+eBM0hdWB+NiHkV9T4D/C8wFvg74GJSt96Xahx7JZLWI70uCyuK1wGuydrxGPBM1q68K4HPZvVmZcc5tOLYnwGuy+qdDWxN+v2ulrXPyhQRvvnWtBswHohubkOzOpcD19c4xqjK+hXlY7LyAdnjodnjqyvqrA+8QXoTX72i/B7gJ908n0gfts4G/lwgjvNJyWXDirJBpLGnE7PHB2T77lpRZ0vgTWBmjbbvlu03ukadG4BHKh4vAF7sel2ysqOBt4G/r/M13S97vGdFnXVJZxhX5p5zPrBGRdmlwMIe/j66nm9g9ppvAfwke112zP0Njc7tOyIr3z57PCx7fEqN3+vjlX8fWfmxwHJgo3b/v/T3m89grAwvA3tXKX8mu58NHCfpPNIg8/0R8VYvnu+Orh8iYomk54BZuWPOA97f9UDSWsBZpDfi9wNrVmxbI7KzrW7sDdwGLKk4U3sFuB/o6kraBVgUEb+viO1xSfc30L4ibotsTCPzS9Ib7M7AI3UcZxfg2YiY1VUQEX+TNBXIX8E2I/c6PQwMlrRmRLzRw/O8VPHz88CxETG7oiyAW3o4xsjsfnI327ch/W5/mjujvpPUzbg96azHSuIEY2V4MyLuq7H9B8B6pK6Vc4EXJF0BfL3BRPNS7vHr3ZStVfF4InA8qdvqD1n90cA5Wb2ldG9j0pnG4VW2dSW7TYFnq2x/ltT27nQN3G9Zo86WFfUqj/uOiFgmaSnVu5Vq2Sx/rMwiUvdTpWqvsYD3ks4ia/kYqWvxeeDJiHg7t31xRLzewzE2Av4WafC/mo2z+5u72b5FD8e3XnKCsZbL3kwuAS7JxhyOBr4FPAVc0aIwDgMmRcQFXQWSDiy474uky4fPr7Ltlex+ITC4yvbBpO6ZqiLiyWwA/pPAv+e3S/oA6ZN3/rkH5+qtAwwgjbfU46/5Y2WGkNrdLH/MnXHlFfkekReAdSWt302S6Yp3LGn8Le8vBZ7DesFXkVlbRcSTEfEdUhfWtllx1yfXtarv1RRrUzHQrTT35Ihcne7iuAPYDngoIu7L3eZmde4FhkjateI53g/0OKEQuAzYS9Inqmz7Vhb3f+bK99G7J0t+ivQm3XUmWfQ1/T2pm+udiwiyZHUg6cKGvuTO7P6z3WyfSzrTG1rl93RfRLzQmjBXXT6DsTKsIWm3KuVPRsTTkq4kfbr8HWm8ZiTwIdJVZZDeGABOUJr3sSwi5jQ5xtuAEyXNy2I5kdS1U6m7OC4GjgHulDSJ9CY2BNgTuDsiriN1yzwA/EzSOFJSOI/q3U95k0jjPL+U9F1gJqlb7TjSYP2/xMpzYJYDN0m6kNTNdSHwy4h4uIe2vEtETJf0G+Anks4knSV8hZSQV5qj0k4RMVfSVcBFkgYDd5Emkn46Io6IiLclnQ78l6T1SWM6rwNbAQdn9Vo9wXTV0u6rDHzrXzdqX0V2TlZnDPBr0hv7MtKlrsfljnM66QqgN4EFFftVu4psVG7fBcB3c2WTgfsqHg8hDYQvIY0vXEC6xPed43cXR1a+OXB1tu9r2XP+CNiuos77SasXLM+OcQLwc2pcRVax7xrAl7PXZjmwmPQGuUeVuguAi7LXfhHwN9KluRvU+5pmZZsAP8yeczlpIHznAq/xSseqEmuROuOB56uUj6DiKrKsbHWyq/9IyeMp4Ae5/fYHfpW9LktIF5l8k4or4Hwr56bsF9Aykj5Iuq5+d1I3w68iYkSB/QaSLoM8mNS1N5V0eeILuXqjSX88HyL90Z0XET9pZhvM+pJszObnEeF5HdantGMMZjvSHIG5pAlURf2U9AnmeNKnoJ1J8wHeobQQ4PXADNKnlpuA67rpyzYzsxK14wxmtcguSZT0c2Djns5gJO0O/IY0+euurGwX0oDkPhFxe1Y2HVgzIj5ese/NwPoR0dAqtGZ9nc9grK9q+RlMrHy9exH7kyat3VVxnHtIlxnuDyDpvaTB4p/m9p0C7N61LpNZfxMRQ51crC/qlMuUhwGPVil/JNsGaY2hNavUe4TUzm1Ki87MzFbSKZcpD2LlWcOQrnLZqqIOVeotzm1/F0ljSROxWHvttXfaYotik3vffvttVlutU/Jzdf2hDeB2lG29x9JQ6Svb9PwZra+2oV5uR3GPPfbY8xGxSbVtnZJgShMRVwFXAQwfPjzuu6/WCicrzJw5kxEjRpQYWfn6QxvA7SidlO7nzq1djz7chjq5HcVJery7bZ2SoheTVl/NG8SKM5Su+3y9QbntZmbWAp2SYB5lxVhLpcqxmfmkBfby9YaRli2v55JoMzPrpU5JMLcAm2bzXADIvmFvq2wbEfEaaf7LYbl9Dwd+GxEvtyhWMzOjDWMw2cJ5B2QP3wesL+nT2eObIy0zPo/0fR7HAUTEbyXdCvxQ0ldIZyQTSes+3V5x+POBmZIuJU3CPCC77Vd6w8zM7F3aMcg/GPhZrqzr8QdIaxytQVpjqNLhpCXef0DFUjGVFSLi7ixZfRP4V9I8maMi4tYmxm+2amnxZGzrP1qeYCJiAelLiWrVGVql7CXg89mt1r43kFtCxsxWNvTMmxred8F3in51jq3KOmUMxszMOowTjJnVdOPkU7lx8qntDsM60Co/0dLMatth0fx2h2AdymcwZmZWCicYMzMrhROMmZmVwgnGzMxK4QRjZmal8FVkZlbTtR/Zt90hWIdygjGzms7e7+R2h2Adyl1kZmZWCicYM6tp+4Xz2H7hvHaHYR3IXWRmVtPUa04DYOi4qW2OxDqNz2DMzKwUTjBmZlYKJxgzMyuFE4yZmZXCCcbMzErhBGNmZqXwZcpmVtOoz13a7hCsQznBmFlND276wXaHYB3KXWRmZlYKJxgzq2nCtElMmDap3WFYB3KCMbOajnpgOkc9ML3dYVgHcoIxM7NSOMGYmVkpnGDMzKwUTjBmZlYKJxgzMyuFJ1qaWU1zhmzd7hCsQznBmFlNB425rN0hWIdyF5mZmZXCCcbMzErhBGNmNS2YOIoFE0e1OwzrQE4wZmZWCicYMzMrhROMmZmVwgnGzMxK4QRjZmalcIIxM7NSeCa/mdV01r4ntTsE61BOMGZW03U77tfuEKxDtbyLTNK2ku6QtEzSM5K+IWn1HvYZLym6uZ1VUW9yN3WGld8yMzOr1NIzGEmDgNuBh4HRwNbARaREd06NXb8PTMuVHQyMA27JlT8KfD5XtqCxiM3syNnpX89nMlavVneRfRFYGzgkIpYAt0laHxgv6YKsbCUR8RTwVGWZpK8Bj0bE7Fz1v0XE70qI3WyV9O3plwNOMFa/VneR7Q9MzyWSKaSks2fRg0jaCNgHuK654ZmZWbO0OsEMI3VhvSMingCWZduKOhRYk+oJZltJSyS9JuluSYUTl5mZNY8ionVPJr0BnBERl+bKnwJ+GBFnFzzOncDAiNgpV34q8DppjGcT4HRgJ2CPiLinm2ONBcYCDBkyZKcpU6YUasvSpUsZMGBAobp9VX9oA7gdjZrz9MuF6p18zMEATPrRDe+U7fC+gVXr+nfRt7SiHSNHjrw/IoZX29ZxlylL2ozUnTYuvy0iLsvVvRl4CDibdFHASiLiKuAqgOHDh8eIESMKxTFz5kyK1u2r+kMbwO1o1JgzbypU7+Ts/qI5K94uFhw9ompd/y76lna3o9VdZIuBah99BmXbivgMIOAnPVWMiGXAzcA/Fg3QzMyao9UJ5lFyYy2StgDWITc2U8MRwN0R8WTB+pHdzMyshVrdRXYLcIak9SLilazscGA5MKunnSUNBXYDvlTkySStDRwI3N9IsGYGQ8dNbXcI1qFafQZzBfAa8AtJe2cD7OOBiysvXZY0T9J/Vtn/COBN4Gf5DZIGSvqVpBMk7SXpcGAGsDkwoYS2mJlZDS09g4mIxZL2Ai4HbgReAi4hJZl8XNWWjzkCuCMinq+y7TXgOdKKAIOBV4HfAntGxH1NaYCZmRXW8qvIIuJh4OM91BnaTfmONfZ5FTikV8GZ2UpunHwqAAeNuayHmmbv1nGXKZtZa+2waH67Q7AO5S8cMzOzUjjBmJlZKZxgzMysFE4wZmZWCicYMzMrha8iM7Oarv3Ivu0OwTqUE4yZ1XT2fif3XMmsCneRmZlZKepKMJKqLd9iZv3Y9gvnsf3Cee0OwzpQvWcwT0u6QNLflxKNmfU5U685janXnNbuMKwD1ZtgrgA+DTwo6feSxkpav4S4zMysw9WVYCJifERsBewDzAUuBv4q6ceS9i4jQDMz60wNDfJHxJ0R8VlgU9JXdn8YmC5pgaTxkjZvZpBmZtZ5ensV2XDgY6SvQV4M/Ao4Hpgn6ZheHtvMzDpY3QlG0paSvi5pPnAHsBlwLLB5RPwLsCVwJXBhUyM1M7OOUtdES0kzgH8GngauBq6OiMcr60TEW5KuBU5tWpRmZtZx6p3J/yxwAHBbRESNerOBDzQclZn1GaM+d2m7Q7AOVW+C+Q/gD9WSi6QBwD9GxF0R8Qbw+Ep7m1nHeXDTD7Y7BOtQ9Y7BzAC27Wbbh7PtZmZmdScY1dg2AFjWi1jMrA+aMG0SE6ZNancY1oF67CKT9DFgREXR8ZL2y1VbCzgQmNO80MysLzjqgemAV1W2+hUZg9mVNJkSIIDDgDdzdV4HHgXOaF5oZmbWyXpMMBFxIdmcFkl/AT4VEbPLDszMzDpbXVeRRYQvPTYzs0KKjMEcANwdEUuyn2uKiJubEpmZmXW0ImcwU4HdgHuyn4PuryYLwF9KZmZmhRLMB4C/VvxsZquQOUO2bncI1qGKDPI/Xu1nM1s1HDTmsnaHYB2qyBjMOvUcMCI82dLMzAp1kS0lja0U5TEYMzMrlGCOpb4EY2b9yIKJowAYOm5qmyOxTlNkDGZyC+IwM7N+prdfmWxmZlZVkUH+e4AxEfGwpHvpobssInZpVnBmZta5iozBPAQsr/jZ4zFmZtajImMwn6/4eUyp0ZiZWb/R8BiMkk0k1foSMjMzW0XVtZoyvLP45TnATtn+b0q6H/hWRNzU5PjMrM3O2vekdodgHaquBCPpBOB7wB3AqcCzwGDgEOB/JH0pIq5sepRm1jbX7Zj/AluzYuo9gzkbuDIivpQrv0LSFcBXAScYMzOrewxmI+CX3Wy7HtiwpwNI2lbSHZKWSXpG0jck1VxeRtJQSVHlNqVK3dGS5kh6VdLDkg4v1DIzq+rI2dM4cva0dodhHajeM5gZwJ7AbVW27QncVWtnSYOA24GHgdHA1sBFpER3ToHn/wrw64rHz+eOvwcp0X0POAU4ALhO0uKIuLXA8c0s59vTLwfcVWb1KzLRctuKh/8OfF/SRsANrBiD+RSwP3B8D4f7IrA2cEhELAFuk7Q+MF7SBVlZLXMj4nc1tn8NuCsiTskez5C0HXAu4ARjZtZCRc5gHuTdkysFnJDd8t9uOY3aqynvD0zPJZIpwETSGdCNBeKpStJ7gZGkM5dKU4CrJQ2MiJcbPb6ZmdWnSIIZ2cTnGwbcWVkQEU9IWpZt6ynBXC1pQ9KZ03XAVyOia5WBrYE1gUdz+zxC6oLbBri3d+GbmVlRRWbyz2ri8w0CXqpSvjjb1p3XgP8gdXMtAUYA40hJZXTFsaly/MW57e8iaSwwFmDIkCHMnDmzVvzvWLp0aeG6fVV/aAO4HY06fYc3G67fXZz+XfQt7W5H3RMtu0haDVgrX17GN1pGxF+BytleMyUtAr4n6SMR8UAvjn0VcBXA8OHDY8SIEYX2mzlzJkXr9lX9oQ3gdjRqzJnF5kWfnN1fNGfF28WCo0dUrevfRd/S7nbUdZlytjzMOEnzgDeAV6rcalkMDKxSPogVZxpF/Ty736ni2FQ5/qDcdjMza4F6z2BOAc4ELgC+BXwTeAs4AngPMKGH/R8ljbW8Q9IWwDqsPHbSk8jdzyclvWFAZbfeMOBt4LE6j29m+JssrXH1TrT8AvB1UoIBuCEizgO2IyWID/Ww/y3AvpLWqyg7nPR1APWO9Xw6u78fICJeI83TOSxX73Dgt76CzMysteo9g/kAMDsi3pL0BrABQES8Lel7wPdJZzjduYJ0FvQLSROBrYDxwMWVly5nXXCzIuK47PF4YD3SJMslwMeAM4BfRMT/Vhz/fNL4zKWkeToHZDfPEDMza7F6z2BeAAZkPz8BfLRi2yDSJMpuRcRiYC/SXJkbgfOAS0hnRZXW4N3zaR4lzZO5GrgZOAq4MLuvPP7dpDObvYHpwCeBozyL36xxN04+lRsnn9ruMKwD1XsG82tgZ9Kb/LWkGfgbAq8DJ5JWWa4pIh4GPt5DnaG5x1NIEyZ7FBE3kM5ezKwJdlg0v90hWIeqN8GMB96X/TyB1EU2hnTmchsrrmg0M7NVXF0JJiLmAnOzn18jfSeMz53NzGwlvZlo+XfAZsAzEfF080IyM7P+oN5BfiT9q6QngceB3wNPSHpKUv5LyMzMbBVW70z+c4HLSfNZDgSGZ/e3AP+ebTczM6u7i+xEYEJEfC1XPi1bG+xE4BtNiczM+oRrP7Jvu0OwDlVvglmb7r+1cha+isys3zl7P/9bW2PqHYO5ATikm22HAl60yMzMgGJfmXxAxcNbgAskDWXlr0zeDvi35odoZu20/cJ5ADy46QfbHIl1miJdZFNZ+auR3wdU65j9EembJs2sn5h6zWmAV1W2+hVJMB8oPQozM+t3inxl8uOtCMTMzPqXumfyS1qDNKC/B7Ah8CLwK9LS+fV9ybeZmfVbdSUYSYOBW4F/ABYAi4DdSfNfHpD0iYh4rtlBmplZ56n3MuWLgY2A3SJiq4jYPSK2AnbNyi9udoBmZtaZ6k0wBwDjIuKeysKIuBc4i7RsjJmZWd1jMO8FXulm2yvAe3oXjpn1NaM+d2m7Q7AOVW+C+R0wTtKdEfG3rkJJ6wLjsu1m1o94gqU1qt4EczowA3hS0q2kQf7BpEmXAkY0NTozM+tYdY3BRMRs4EPAVcAmwD6kBHMF8KGIeKDpEZpZW02YNokJ0ya1OwzrQIXPYCStCewC/CUiziwvJDPrS456YDrgVZWtfvWcwbwF3AkMKykWMzPrRwonmIh4G/gTsGl54ZiZWX9R7zyYrwLnStqhjGDMzKz/qPcqsnNIM/ZnS3qadBVZVFaIiF2aFJuZmXWwehPMg9nNzMyspkIJRtLapGViHgQWArdHxKIyAzOzvmHOkK3bHYJ1qCJfmbwVcDswtKJ4iaTPRMStZQVmZn3DQWMua3cI1qGKDPJfALwN/DOwDrAd8EfgyhLjMjOzDlckwewOnBMRv46IVyPiEeAE4P2SNis3PDMz61RFEsxmwJ9zZfNJa495ToxZP7dg4igWTBzV7jCsAxWdBxM9VzEzM1uh6GXK0yW9WaX8jnx5RAzufVhmZtbpiiSY80qPwszM+p0eE0xEOMGYmVnd6l2LzMzMrBAnGDMzK0W9a5GZ2SrmrH1PancI1qGcYMysput23K/dIViHcheZmZmVwgnGzGo6cvY0jpw9rd1hWAdqeYKRtK2kOyQtk/SMpG9IWr2HfXaWdLWkedl+cyV9XdJauXrjJUWVm8/xzRr07emX8+3pl7c7DOtALR2DkTSItPT/w8BoYGvgIlKiO6fGrodndScCfwL+ATg/uz80V/dlIJ9QHult7GZmVp9WD/J/EVgbOCQilgC3SVofGC/pgqysmu9ExPMVj2dKehW4UtKWEfF4xbY3I+J35YRvZmZFtbqLbH9gei6RTCElnT272ymXXLr8MbvfvHnhmZlZs7Q6wQwDHq0siIgngGXZtnrsTvoitPm58g0kPS/pDUl/lHRIw9GamVnDFNG6lfglvQGcERGX5sqfAn4YEWcXPM6mwP8CN0fEmIryY4DBpLOb9UhfjHYAcGhE/KKbY40FxgIMGTJkpylTphRqy9KlSxkwYEChun1Vf2gDuB2NmvP0y4XqnXzMwQBM+tEN75Tt8L6BVev6d9G3tKIdI0eOvD8ihlfb1nETLSW9B/gpsBT4cuW2iPhRru6NwG+Ac4GqCSYirgKuAhg+fHiMGDGiUBwzZ86kaN2+qj+0AdyORo0586ZC9U7O7i+as+LtYsHRI6rW9e+ib2l3O1qdYBYD1T76DMq21SRJwA+B7YB/ioia+0RESPoFMFHS6hHxVgMxm63Sho6b2u4QrEO1OsE8Sm6sRdIWwDrkxma6cSnp8uZ9IqJIfUjfxulv5DQza7FWD/LfAuwrab2KssOB5cCsWjtKOgs4CTgmIu4u8mTZGc+hwAM+ezEza61Wn8FcAZwC/ELSRGArYDxwceWly5LmAbMi4rjs8VHABGAy8LSk3SqOOT8insvqzQKuJ50NrQt8AdgVOLjcZpn1XzdOPhWAg8Zc1uZIrNO0NMFExGJJewGXAzcCLwGXkJJMPq7K5WM+kd2PyW6VPk9KPADzgNOAzUiXMP8BODAibmlG/Garoh0W5WcCmBXT8qvIIuJh4OM91BmaezyGlRNLtf2O60VoZmbWRF5N2czMSuEEY2ZmpXCCMTOzUjjBmJlZKTpuqRgza61rP7Jvu0OwDuUEY2Y1nb3fyT1XMqvCXWRmZlYKJxgzq2n7hfPYfuG8dodhHchdZGZW09RrTgO8qrLVz2cwZmZWCicYMzMrhROMmZmVwgnGzMxK4QRjZmalcIIxM7NS+DJlM6tp1OcubXcI1qGcYMyspgc3/WC7Q7AO5S4yMzMrhROMmdU0YdokJkyb1O4wrAM5wZhZTUc9MJ2jHpje7jCsAznBmJlZKZxgzMysFE4wZmZWCicYMzMrhROMmZmVwhMtzaymOUO2bncI1qGcYMyspoPGXNbuEKxDuYvMzMxK4QRjZmalcIIxs5oWTBzFgomj2h2GdSAnGDMzK4UTjJmZlcIJxszMSuEEY2ZmpXCCMTOzUjjBmJlZKTyT38xqOmvfk9odgnUoJxgzq+m6HfdrdwjWodxFZmZmpXCCMbOajpw9jSNnT2t3GNaB3EVmZjV9e/rlgLvKrH4+gzEzs1K0PMFI2lbSHZKWSXpG0jckrV5gv4GSrpa0WNLLkn4saaMq9UZLmiPpVUkPSzq8nJaYmVktLU0wkgYBtwMBjAa+AZwOnFdg958CI4DjgTHAzsANuePvAVwPzAD2B24CrpP0iaY0wMzMCmv1GMwXgbWBQyJiCXCbpPWB8ZIuyMpWIml34BPAnhFxV1b2NPB7SXtHxO1Z1a8Bd0XEKdnjGZK2A84Fbi2vWWZWy9Azb2p43wXfObCJkVgrtTrB7A9MzyWSKcBEYE/gxhr7LepKLgARcY+kv2Tbbpf0XmAkcEpu3ynA1ZIGRsTLTWqHWdPNefplxjTwRuw3YOurWp1ghgF3VhZExBOSlmXbuksww4BHq5Q/km0D2BpYs0q9R0hdgdsA9zYWtq1q/Im78/h31ve0OsEMAl6qUr4429bIfltV1KFKvcW57e8iaSwwNnu4VNLcGnFU2hh4vmDdvqo/tAH6WDs0seFdG2pHL56v2PG7fqj4Vssaz9n030XZ7evmOfrU31QvtKIdW3a3YZWfBxMRVwFX1bufpPsiYngJIbVMf2gDuB19SX9oA7gdzdLqy5QXAwOrlA9ixZlGo/t13efrDcptNzOzFmh1gnmUFWMmAEjaAliH6mMs3e6XqRybmQ+8UaXeMOBt4LEG4jUzswa1OsHcAuwrab2KssOB5cCsHvbbNJvnAoCk4aTxl1sAIuI10vyXw3L7Hg78toQryOruVuuD+kMbwO3oS/pDG8DtaApFROueLE20fBh4kHRp8lbAxcClEXFORb15wKyIOK6ibDrwIeArpDOSicCzEfHPFXX2AGYCl5MmYR6Q1d8vIjwPxsyshVp6BhMRi4G9gNVJlySfB1wCfD1XdY2sTqXDSWc5PwB+CNwPfCp3/LuBTwN7A9OBT2bJikEAAAZPSURBVAJHObmYmbVeS89gzMxs1eHVlAuS9AVJf8oW0bxf0l517v9RSW9Jauu19Y20Q9IJkm6TtChbaPTXrVjfreyFUVulkXZI2jlrw7xsv7mSvi5prVbFXSWmhn4fFfuvJuk+SSFpVM97NF9v2iDpEEn3Slou6QVJ0yStW3bM3cTS6P/GcEm3Snoxu90uadfSAo0I33q4AUcCb5HWOhtJ6qJbDmxfcH8BvwYWAs93WjuAJ0iDhQcD+wDXkMbBPllirIOAZ0iLo+5DWsfub8A3C+w7HfgLcCipG/Ux4Fdtes0bagfwXeAu4AukRV5PAV4Gru+kduSOMTb7HwhgVCe1gbTI7qukBXpHZH9Xk4CBndIOYAvSRPQ7gQOz20xgCbBlKbG2+sXpxBswF/hBxePVgDnAjwru/y/APGBCmxNMQ+0ANq5S9htgRomxnkWau7R+Rdm/Acsqy6rst3v2BvaxirJdsrK92/CaN9qOaq/52KwdW3ZKOyrqDgKeA45rY4Jp+HcBvAJ8odUxN7kdXyR9wBxYUTYoK/vXMmJ1F1kPJG1FWsfsp11lEfE28DPSQps97b8e6Yq3rwCvlxRmj3rTjoio1q33R2DzZsaY093CqGuTFkattd9KC6OSzmh6/H2VoKF21HjNodzXvTuN/j66nE86i7+jhNiKarQNn8nurykrsDo12o41gTdJZztdlmZlqrpHLznB9Kxr4ma1RTQ3lLRJD/ufCzwSETf0UK9svW1H3u6UO3l1pQVOI+IJ0qe0apNuu90vU7kwais12o5qdid1Tc5vTmh1abgdkv4BOJb0IaudGm3DrqSz/+MkPSXpDUm/l/R/ygu1pkbbcX1W5yJJgyUNJl3Fu5j0QbPpnGB61tAimgCSPgycCJxWQlz1argdeZKOBT5KmsNUljIWRi3cxiZqSjySNgXOAf4rIp5tUmz16E07JgGXR8S8pkdVn0bbsCnwYdLrPw44iHQWME3SkGYHWUBD7YiIZ0hjr4cCi7LbIcC+EfFcCXGumotdShoIbNZTvYiotXxNEZcBkyNiTi+PU1UL21H5nDuR3jAui4gZzTqudU/Se0hdm0uBL7c5nLpIOoL05nxQu2PpBQEDgMMiYhqApN8AjwMnkS6a6fMkbUY6U7mfdNECpA/AN0n6P9lZUFOtkgmGtJzM/ytQT7x7Ec3KTw01F9GUtD/wT8BJkjbIitdKm7QBsDzS8ja9UXo73nWQNI5zE6kf/fTiYTakNwujVuvu62m/sjTaDiD9sZCu9tsO+KdIk5Xboe52SFoTuJA0Brla9ne/frZ5XUnrRcQrZQTbjd78TQXpiisAImKJpPuBbZsZYEGNtuMM0jjMpyPiDQBJdwJ/InVf5r+ssddWyS6yiPh+RKinW1a969N/tUU0X6xxavlh0qeeP5F+6YtJp9cbZj+f0SHtACDrr51O+tR2RES81dv4e1Dmwqit1Gg7ulwKjAZGN/NMtAGNtGNd4O9IXald/wMPZNumsOKihVZp9HfxCOlDWn4gXKQxsVZrtB3DgIe6kgtARLwOPET6wsamWyUTTD0i4s+kwex3FtGUtFr2+JYau/6c1N9ZebuGdM35SOC/Sgq5ql60A0kDgJuzh6MiYllZcVYobWHUFmu0HUg6i9QFc0ykZZDaqZF2LGXl/4Ejs21nA0eXE2q3Gv1dTM3uR3YVZN3TO7EiYbZSo+14HNg+63IFQOmr5rcHFpQQp+fBFLmxYoLiOaQ/ssnkJiiSLg98E9izxnHG0zcmWtbVDuBW0iXWRwG7Vd5KjHUQ8FfgNtLacmNJb1jfzNWbB/xnrmw68GfSAObBpCuA2jnRsu52ZK91AFfnX3Ngk05pR5XjDKW9Ey0b/Zu6Idv3c6QJirNI83oGdUo7SAnxDVI394HAKFKyegP4SCmxtvrF6dQbaUb1POA14A/AXrntI7J/nBE1jjGeNiaYRtuRPa56KznWbUmzjpdn/1DnA6vn6iwgXUhRWbZB9sb8EumM8VqqTFxs4WtedztIyb+7131Mp7SjyjHalmB6+Tc1APi/wAvZvrcDO3TS31RWthdphYgXs9usWu9Zvb15sUszMyuFx2DMzKwUTjBmZlYKJxgzMyuFE4yZmZXCCcbMzErhBGNmZqVwgjEzs1I4wZiZWSn+P1iuTBNlR9GZAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "# plot estimated values for \"a\"\n",
    "plt.bar(result['values'], result['probabilities'], width=0.5/len(result['probabilities']))\n",
    "plt.xticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.yticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.title('\"a\" Value', size=15)\n",
    "plt.ylabel('Probability', size=15)\n",
    "plt.ylim((0,1))\n",
    "plt.grid()\n",
    "plt.show()\n",
    "\n",
    "# plot estimated values for option price (after re-scaling and reversing the c_approx-transformation)\n",
    "plt.bar(result['mapped_values'], result['probabilities'], width=1/len(result['probabilities']))\n",
    "plt.plot([exact_value, exact_value], [0,1], 'r--', linewidth=2)\n",
    "plt.xticks(size=15)\n",
    "plt.yticks([0, 0.25, 0.5, 0.75, 1], size=15)\n",
    "plt.title('Estimated Option Price', size=15)\n",
    "plt.ylabel('Probability', size=15)\n",
    "plt.ylim((0,1))\n",
    "plt.grid()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
